BEI and equity paper pre-analyses

Data extraction and pre-analyses for the paper looking at BEI and equity. BEI comprise bike lanes and canopy changes while equity is measured through Pampalon deprivation index. (Partially adapted from original work on BEI bike lanes – see bike_lane_stats.R and ReadMe.md.)

Paper is available here.

General processing steps:

• Get CT 2016 boundaries (from Cancensus API)
• Get BEI interventions (canopy / bike lanes) for the selected years (2011 / 2015 / 2016 / 2019)
• Compute changes between years
  • for the original CT boundaries
  • for buffered CTs (250 / 500 / 750m)

In a second phase, these BEI changes will be linked to Pampalon index for 2016 (and 2011 ?)

1 Built Environment Intervention Extraction

1.1 Bike lane changes

We use data categorized by Philippe Apparicio’s team who manually identified bike lanes for each census year since 1991. For this study, we limit ourselves to 2016 and 2011 census years.

On top of the original CT boundaries, three levels of buffer have been applied to the CT – 250m, 500m & 750m. Then the same series of processing steps (see above) have been applied to the buffers.

UPDATE 2021-12-02 Following discussion with @Yan, add normalized bike line changes:

• by CT/buffer area
• by street length within CT/buffer

# Bike lanes, from Ph. Apparicio
reseau <- st_read(dsn="data/ReseauCyclableFinal.gdb", layer = "Reseau") # Already in NAD83 / MTM zone 8

## Reading layer `Reseau' from data source 
##   `/Users/benoit/WORKSPACE/gentrification_BEI_equity/data/ReseauCyclableFinal.gdb' 
##   using driver `OpenFileGDB'
## Simple feature collection with 82166 features and 72 fields
## Geometry type: GEOMETRY
## Dimension:     XYZ, XYZM
## Bounding box:  xmin: 266985.5 ymin: 5029251 xmax: 320986.1 ymax: 5062652
## z_range:       zmin: 0 zmax: 43
## m_range:       mmin: 0 mmax: 43
## Projected CRS: NAD83 / MTM zone 8

bike_lane <- reseau %>%
  filter(An2016 == 1 | An2011 == 1) %>%
  select(IdRte, ClsRte, Zone, starts_with("An"), starts_with("Typo_")) %>%
  st_cast("MULTILINESTRING") # Get rid of a few MULTICURVE geometries

# CT boundaries for Montreal
CT16 <- get_census(dataset='CA16', regions=list(CMA='24462'), level='CT', geo_format = "sf") %>%
  filter(Type == "CT") %>%
  mutate(interact_aoi = CD_UID %in% c(2466, 2465, 2458)) %>% # Flag Montréal island, Laval and the South shore (Longueuil, St-Lambert, Brossard)
  st_transform(st_crs(bike_lane))

## Reading geo data from local cache.

CT11 <- get_census(dataset='CA11', regions=list(CMA='24462'), level='CT', geo_format = "sf") %>%
  filter(Type == "CT") %>%
  st_transform(st_crs(bike_lane))

## Reading geo data from local cache.

# Load gentrified CTs, 5 year span (from repo gentrification_metrics)
ding <- list()
ding[["2016"]] <- st_read("~/WORKSPACE/gentrification_metrics/R/data/gentrified_5years.gpkg", "gentrified_ding_16", quiet=TRUE) %>%
  filter(cma_uid_16 == "24462") %>%
  st_transform(st_crs(bike_lane))
ding[["2011"]] <- st_read("~/WORKSPACE/gentrification_metrics/R/data/gentrified_5years.gpkg", "gentrified_ding_11", quiet=TRUE) %>%
  filter(cma_uid_11 == "24462") %>%
  st_transform(st_crs(bike_lane))
ding[["2006"]] <- st_read("~/WORKSPACE/gentrification_metrics/R/data/gentrified_5years.gpkg", "gentrified_ding_06", quiet=TRUE) %>%
  filter(cma_uid_06 == "24462") %>%
  st_transform(st_crs(bike_lane))

compute_bikelane_by_area <- function(sf_areas, year_fld, typo_fld) {
  # Compute length of bike lanes within each area.
  # --
  # Parameters:
  #   - sf_areas: sf class object defining the areas of interest, must have a GeoUID field
  #   - year_fld: field name specifying the year of interest, e.g. An2016
  #   - typo_fld: field name specifying the typology for the year of interest, e.g. Typo_2016
  
  year_fld <- enquo(year_fld)
  typo_fld <- enquo(typo_fld)
  
  # Compute intersection of bike lanes with areas
  bk <- bike_lane %>%
    filter(!!year_fld == 1) %>%
    st_intersection(sf_areas) %>%
    mutate(bike_lane_length = st_length(.)) %>%
    as.data.frame() %>%
    group_by(GeoUID, !!typo_fld) %>%
    summarise(bike_lane_length = sum(bike_lane_length)) %>%
    ungroup() %>%
    pivot_wider(names_from = !!typo_fld, names_prefix = "Bike_class", names_sort = TRUE,
                values_from = bike_lane_length, values_fill = units::set_units(0, m)) %>%
    mutate(Bike_lane_total = units::set_units(rowSums(select(., starts_with("Bike_class"))), m))
  
  # Merge back into original sf_areas
  bk <- sf_areas %>%
    left_join(bk)
  
  # Replace NA by 0, which occur in Bike_class length
  bk[is.na(bk)] <- 0
  
  return(bk)
}

compute_streetlength_by_area <- function(sf_areas) {
  # Compute length of streets within each area.
  # --
  # Parameters:
  #   - sf_areas: sf class object defining the areas of interest, must have a GeoUID field

  # Compute intersection of streets with areas
  bk <- reseau  %>%
    st_cast("MULTILINESTRING") %>% # Get rid of a few MULTICURVE geometries
    st_intersection(sf_areas) %>%
    mutate(street_length = st_length(.)) %>%
    as.data.frame() %>%
    group_by(GeoUID) %>%
    summarise(street_length = sum(street_length)) %>%
    ungroup()
  
  # Merge back into original sf_areas
  bk <- sf_areas %>%
    mutate(shape_area_km2 = units::set_units(st_area(.), 'km^2')) %>%
    left_join(bk)
  
  # Replace NA by 0
  bk[is.na(bk)] <- 0
  
  return(bk)
}


# Compute year 2016 and year 2011 bike lanes within 2016 CTs
# NB: contrary to the original work, we keep the same area of reference, i.e. 2016
bike_lane_by_CT16 <- compute_bikelane_by_area(CT16, An2016, Typo_2016)

## Warning: attribute variables are assumed to be spatially constant throughout all
## geometries

bike_lane_by_CT11 <- compute_bikelane_by_area(CT16, An2011, Typo_2011)

## Warning: attribute variables are assumed to be spatially constant throughout all
## geometries

# Compute buffers, with 3 radii and for each census year
radii <- c(250, 500, 750)

buf_CT16 <- lapply(radii, st_buffer, x=CT16)
names(buf_CT16) <- lapply(radii, function(b) {paste0("buf", b)})

# Compute bike length for each buffer/census year
buf_CT16_w_bike_length <- lapply(buf_CT16, compute_bikelane_by_area, year_fld=An2016, typo_fld=Typo_2016)

## Warning: attribute variables are assumed to be spatially constant throughout all
## geometries

## Warning: attribute variables are assumed to be spatially constant throughout all
## geometries

## Warning: attribute variables are assumed to be spatially constant throughout all
## geometries

names(buf_CT16_w_bike_length) <- lapply(radii, function(b) {paste0("buf", b)})
buf_CT16_w_bike_length$original <- bike_lane_by_CT16

buf_CT11_w_bike_length <- lapply(buf_CT16, compute_bikelane_by_area, year_fld=An2011, typo_fld=Typo_2011)

## Warning: attribute variables are assumed to be spatially constant throughout all
## geometries

## Warning: attribute variables are assumed to be spatially constant throughout all
## geometries

## Warning: attribute variables are assumed to be spatially constant throughout all
## geometries

names(buf_CT11_w_bike_length) <- lapply(radii, function(b) {paste0("buf", b)})
buf_CT11_w_bike_length$original <- bike_lane_by_CT11

# Compute total street length within CT/buffer
street_length_by_CT16 <- compute_streetlength_by_area(CT16)

## Warning: attribute variables are assumed to be spatially constant throughout all
## geometries

buf_CT16_street_length <- lapply(buf_CT16, compute_streetlength_by_area)

## Warning: attribute variables are assumed to be spatially constant throughout all
## geometries

## Warning: attribute variables are assumed to be spatially constant throughout all
## geometries

## Warning: attribute variables are assumed to be spatially constant throughout all
## geometries

names(buf_CT16_street_length) <- lapply(radii, function(b) {paste0("buf", b)})
buf_CT16_street_length$original <- street_length_by_CT16

# Reorganize data to have all data in one dataframe
bike_lane_changes <- CT16 %>%
  left_join(select(as.data.frame(buf_CT16_street_length$original), GeoUID, street_length, shape_area_km2), by="GeoUID") %>%
  left_join(select(as.data.frame(buf_CT16_street_length$buf250), GeoUID, street_length, shape_area_km2), by="GeoUID", suffix=c(".ct", ".b250")) %>%
  left_join(select(as.data.frame(buf_CT16_street_length$buf500), GeoUID, street_length, shape_area_km2), by="GeoUID") %>%
  left_join(select(as.data.frame(buf_CT16_street_length$buf750), GeoUID, street_length, shape_area_km2), by="GeoUID", suffix=c(".b500", ".b750")) %>%
  left_join(select(as.data.frame(buf_CT16_w_bike_length$original), GeoUID, starts_with("Bike_")), by="GeoUID") %>%
  left_join(select(as.data.frame(buf_CT11_w_bike_length$original), GeoUID, starts_with("Bike_")), by="GeoUID", suffix=c(".2016ct", ".2011ct")) %>%
  left_join(select(as.data.frame(buf_CT16_w_bike_length$buf250), GeoUID, starts_with("Bike_")), by="GeoUID") %>%
  left_join(select(as.data.frame(buf_CT11_w_bike_length$buf250), GeoUID, starts_with("Bike_")), by="GeoUID", suffix=c(".2016b250", ".2011b250")) %>%
  left_join(select(as.data.frame(buf_CT16_w_bike_length$buf500), GeoUID, starts_with("Bike_")), by="GeoUID") %>%
  left_join(select(as.data.frame(buf_CT11_w_bike_length$buf500), GeoUID, starts_with("Bike_")), by="GeoUID", suffix=c(".2016b500", ".2011b500")) %>%
  left_join(select(as.data.frame(buf_CT16_w_bike_length$buf750), GeoUID, starts_with("Bike_")), by="GeoUID") %>%
  left_join(select(as.data.frame(buf_CT11_w_bike_length$buf750), GeoUID, starts_with("Bike_")), by="GeoUID", suffix=c(".2016b750", ".2011b750"))

# Compute change between 2011 and 2016 (only for total bike lane length)
bike_lane_changes <- bike_lane_changes %>%
  mutate(Bike_lane_diff.2011.2016ct = Bike_lane_total.2016ct - Bike_lane_total.2011ct,
         Bike_lane_diff.2011.2016b250 = Bike_lane_total.2016b250 - Bike_lane_total.2011b250,
         Bike_lane_diff.2011.2016b500 = Bike_lane_total.2016b500 - Bike_lane_total.2011b500,
         Bike_lane_diff.2011.2016b750 = Bike_lane_total.2016b750 - Bike_lane_total.2011b750)

# Normalize bike lane change by (i) street length and (ii) area
bike_lane_changes <- bike_lane_changes %>%
  mutate(Bike_lane_diff.by.street.2011.2016ct = Bike_lane_diff.2011.2016ct / street_length.ct,
         Bike_lane_diff.by.street.2011.2016b250 = Bike_lane_diff.2011.2016b250 / street_length.b250,
         Bike_lane_diff.by.street.2011.2016b500 = Bike_lane_diff.2011.2016b500 / street_length.b500,
         Bike_lane_diff.by.street.2011.2016b750 = Bike_lane_diff.2011.2016b750 / street_length.b750,
         Bike_lane_diff.by.area.2011.2016ct = Bike_lane_diff.2011.2016ct / shape_area_km2.ct,
         Bike_lane_diff.by.area.2011.2016b250 = Bike_lane_diff.2011.2016b250 / shape_area_km2.b250,
         Bike_lane_diff.by.area.2011.2016b500 = Bike_lane_diff.2011.2016b500 / shape_area_km2.b500,
         Bike_lane_diff.by.area.2011.2016b750 = Bike_lane_diff.2011.2016b750 / shape_area_km2.b750)

# Save results
st_write(bike_lane_changes, dsn = "data/bike_length_changes.gpkg", delete_layer = TRUE)

## Deleting layer `bike_length_changes' using driver `GPKG'
## Writing layer `bike_length_changes' to data source 
##   `data/bike_length_changes.gpkg' using driver `GPKG'
## Writing 970 features with 96 fields and geometry type Multi Polygon.

# Clean up
rm(buf_CT16)
rm(bike_lane_by_CT11, bike_lane_by_CT16)
rm(buf_CT11_w_bike_length, buf_CT16_w_bike_length)
rm(street_length_by_CT16, buf_CT16_street_length)

Check output for one specific dataset (Census tracts 2016, no buffer)

ggplot() +
   geom_sf(data=filter(bike_lane_changes, interact_aoi), mapping = aes(fill=as.numeric(Bike_lane_total.2016ct)), lwd=0) +
  scale_fill_continuous(name = "Total length (m)")+ 
  labs(title = "Length of bike lanes within 2016 CTs", subtitle = "(INTERACT study area || for control only)")

1.1.1 Basic stats on each layer

1.1.1.1 Census Tracts

# CT level
ggplot(bike_lane_changes) +
  geom_histogram(aes(as.numeric(Bike_lane_diff.2011.2016ct)), binwidth = 250) + 
  xlab("Difference of bike lane length between 2016 and 2011 | CT level")

summary(bike_lane_changes$Bike_lane_diff.2011.2016ct)

##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
## -1377.7     0.0     0.0   262.1   189.7 14463.8

1.1.1.2 Buffers 250m

# buf250 level
ggplot(bike_lane_changes) +
  geom_histogram(aes(as.numeric(Bike_lane_diff.2011.2016b250)), binwidth = 250) + 
  xlab("Difference of bike lane length between 2016 and 2011 | buf 250m")

summary(bike_lane_changes$Bike_lane_diff.2011.2016b250)

##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
## -3452.5     0.0     0.0   646.7   933.8 20786.8

1.1.1.3 Buffers 500m

# buf500 level
ggplot(bike_lane_changes) +
  geom_histogram(aes(as.numeric(Bike_lane_diff.2011.2016b500)), binwidth = 250) + 
  xlab("Difference of bike lane length between 2016 and 2011 | buf 500m")

summary(bike_lane_changes$Bike_lane_diff.2011.2016b250)

##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
## -3452.5     0.0     0.0   646.7   933.8 20786.8

1.1.1.4 Buffers 750m

# buf750 level
ggplot(bike_lane_changes) +
  geom_histogram(aes(as.numeric(Bike_lane_diff.2011.2016b750)), binwidth = 250) + 
  xlab("Difference of bike lane length between 2016 and 2011 | buf 750m")

summary(bike_lane_changes$Bike_lane_diff.2011.2016b750)

##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
## -3697.5     0.0   578.6  1842.2  2928.6 29325.1

1.2 Canopy changes

Canopy changes is based on data produced by CMM, using multispectral aerial imagery and lidar. In order to sync the observations with the census years, we focus on 2011 and 2019 with one extra observation point in 2015.

The processing steps are similar to the ones for the bike lanes:

• Import the raster for each of the 3 years
• Compute proportion of canopy within each area level (CT and buffers) for the 3 years

# Codes du raster "espace vert"
# 0. No data (hors CMM)
# 1. NDVI < 0,3 et MNH < 3,0m = Minéral bas (route, stationnement, etc.)
# 2. NDVI < 0,3 et MNH ≥ 3,0m = Minéral haut (constructions)
# 3. NDVI ≥ 0,3 et MNH < 3,0m = Végétal bas (culture, gazon, etc.)
# 4. NDVI ≥ 0,3 et MNH ≥ 3,0m = Végétal haut (canopée)
# 5. Aquatique

# Load rasters into pg database for further processing
system("psql -d xgentrif_bei -c 'CREATE EXTENSION IF NOT EXISTS postgis'")
system("psql -d xgentrif_bei -c 'CREATE EXTENSION IF NOT EXISTS postgis_raster'")
if (nrow(dbGetQuery(con_bei, "SELECT 1 test WHERE to_regclass('canopee2019') IS NOT NULL;")) == 0) {
  system("raster2pgsql -s 32188 -I -C -M data/canopy/2019/*.tif -F -t 1000x1000 canopee2019 | psql -d xgentrif_bei", intern = TRUE)
} else { message("PG Raster 'canopee2019' already imported") }

## PG Raster 'canopee2019' already imported

if (nrow(dbGetQuery(con_bei, "SELECT 1 test WHERE to_regclass('canopee2015') IS NOT NULL;")) == 0) {
  system("raster2pgsql -s 32188 -I -C -M data/canopy/2015/*.tif -F -t 1000x1000 canopee2015 | psql -d xgentrif_bei", intern = TRUE)
} else { message("PG Raster 'canopee2015' already imported") }

## PG Raster 'canopee2015' already imported

if (nrow(dbGetQuery(con_bei, "SELECT 1 test WHERE to_regclass('canopee2011') IS NOT NULL;")) == 0) {
  system("raster2pgsql -s 32188 -I -C -M data/canopy/2011/*.tif -F -t 1000x1000 canopee2011 | psql -d xgentrif_bei", intern = TRUE)
} else { message("PG Raster 'canopee2011' already imported") }

## PG Raster 'canopee2011' already imported

# Resample to 10m as the original rasters have a 1m resolution, which is too high to allow for a swift processing
if (nrow(dbGetQuery(con_bei, "SELECT 1 test WHERE to_regclass('canopee2019_10m') IS NOT NULL;")) == 0) {
  system("gdal_translate -of GTiff PG:\"host=localhost dbname=xgentrif_bei table=canopee2019 mode=2\" -r mode -tr 10 10 data/canopy/canopee2019_10m.tif")
  system("raster2pgsql -s 32188 -I -C -M data/canopy/canopee2019_10m.tif -F -t 100x100 canopee2019_10m | psql -d xgentrif_bei")
} else { message("PG Raster 'canopee2019_10m' already imported") }

## PG Raster 'canopee2019_10m' already imported

if (nrow(dbGetQuery(con_bei, "SELECT 1 test WHERE to_regclass('canopee2015_10m') IS NOT NULL;")) == 0) {
  system("gdal_translate -of GTiff PG:\"host=localhost dbname=xgentrif_bei table=canopee2015 mode=2\" -r mode -tr 10 10 data/canopy/canopee2015_10m.tif")
  system("raster2pgsql -s 32188 -I -C -M data/canopy/canopee2015_10m.tif -F -t 100x100 canopee2015_10m | psql -d xgentrif_bei")
} else { message("PG Raster 'canopee2015_10m' already imported") }

## PG Raster 'canopee2015_10m' already imported

if (nrow(dbGetQuery(con_bei, "SELECT 1 test WHERE to_regclass('canopee2011_10m') IS NOT NULL;")) == 0) {
  system("gdal_translate -of GTiff PG:\"host=localhost dbname=xgentrif_bei table=canopee2011 mode=2\" -r mode -tr 10 10 data/canopy/canopee2011_10m.tif")
  system("raster2pgsql -s 32188 -I -C -M data/canopy/canopee2011_10m.tif -F -t 100x100 canopee2011_10m | psql -d xgentrif_bei")
} else { message("PG Raster 'canopee2011_10m' already imported") }

## PG Raster 'canopee2011_10m' already imported

# Push CT16 to pg
if (nrow(dbGetQuery(con_bei, "SELECT 1 test WHERE to_regclass('ct16') IS NOT NULL;")) == 0) {
  CT16 %>%
    st_transform(crs = 32188) %>%
    st_write(con_bei, "ct16",
             layer_options = c("OVERWRITE=yes", "LAUNDER=true", "SPATIAL_INDEX=gist", "GEOMETRY_NAME=geom"))
  system("psql -d xgentrif_bei -c 'CREATE INDEX ON  ct16 USING gist (geometry)'")
} else { message("PG Layer CT16 already imported") }

## PG Layer CT16 already imported

1.2.1 Extract % of green spaces at various scales

1.2.1.1 Census Tracts

WITH cnt19 AS (
    SELECT "GeoUID"
        ,(pvc).value, SUM((pvc).count) As total
    FROM (SELECT "GeoUID"
            ,ST_ValueCount(ST_Clip(rast, geometry)) As pvc
        FROM canopee2019_10m
        JOIN ct16 ON ST_Intersects(geometry, rast)
    ) As foo
    GROUP BY "GeoUID", (pvc).value
),
canopee19 AS (
    SELECT "GeoUID"
        ,round(100. * sum(total) FILTER (WHERE value in (3, 4)) / sum(total), 2) AS pct_esp_vert_2019
    FROM cnt19
    WHERE value > 0 -- discard no data, including postgis raster no data
    GROUP BY "GeoUID"
),
cnt15 AS (
    SELECT "GeoUID"
        ,(pvc).value, SUM((pvc).count) As total
    FROM (SELECT "GeoUID"
            ,ST_ValueCount(ST_Clip(rast, geometry)) As pvc
        FROM canopee2015_10m
        JOIN ct16 ON ST_Intersects(geometry, rast)
    ) As foo
    GROUP BY "GeoUID", (pvc).value
),
canopee15 AS (
    SELECT "GeoUID"
        ,round(100. * sum(total) FILTER (WHERE value in (3, 4)) / sum(total), 2) AS pct_esp_vert_2015
    FROM cnt15
    WHERE value > 0 -- discard no data, including postgis raster no data
    GROUP BY "GeoUID"
),
cnt11 AS (
    SELECT "GeoUID"
        ,(pvc).value, SUM((pvc).count) As total
    FROM (SELECT "GeoUID"
            ,ST_ValueCount(ST_Clip(rast, geometry)) As pvc
        FROM canopee2011_10m
        JOIN ct16 ON ST_Intersects(geometry, rast)
    ) As foo
    GROUP BY "GeoUID", (pvc).value
),
canopee11 AS (
    SELECT "GeoUID"
        ,round(100. * sum(total) FILTER (WHERE value in (3, 4)) / sum(total), 2) AS pct_esp_vert_2011
    FROM cnt11
    WHERE value > 0 -- discard no data, including postgis raster no data
    GROUP BY "GeoUID"
)
SELECT "GeoUID"
    ,st_area(geometry) ct_area_m2
    ,coalesce(pct_esp_vert_2011, 0) pct_esp_vert_2011
    ,coalesce(pct_esp_vert_2015, 0) pct_esp_vert_2015
    ,coalesce(pct_esp_vert_2019, 0) pct_esp_vert_2019
FROM ct16
FULL JOIN canopee19 USING ("GeoUID")
FULL JOIN canopee15 USING ("GeoUID")
FULL JOIN canopee11 USING ("GeoUID");

1.2.1.2 Buffers 250m

WITH ct16 AS (
    select "GeoUID", ST_Buffer(geometry, 250) geometry
    from ct16
),
cnt19 AS (
    SELECT "GeoUID"
        ,(pvc).value, SUM((pvc).count) As total
    FROM (SELECT "GeoUID"
            ,ST_ValueCount(ST_Clip(rast, geometry)) As pvc
        FROM canopee2019_10m
        JOIN ct16 ON ST_Intersects(geometry, rast)
    ) As foo
    GROUP BY "GeoUID", (pvc).value
),
canopee19 AS (
    SELECT "GeoUID"
        ,round(100. * sum(total) FILTER (WHERE value in (3, 4)) / sum(total), 2) AS pct_esp_vert_2019
    FROM cnt19
    WHERE value > 0 -- discard no data, including postgis raster no data
    GROUP BY "GeoUID"
),
cnt15 AS (
    SELECT "GeoUID"
        ,(pvc).value, SUM((pvc).count) As total
    FROM (SELECT "GeoUID"
            ,ST_ValueCount(ST_Clip(rast, geometry)) As pvc
        FROM canopee2015_10m
        JOIN ct16 ON ST_Intersects(geometry, rast)
    ) As foo
    GROUP BY "GeoUID", (pvc).value
),
canopee15 AS (
    SELECT "GeoUID"
        ,round(100. * sum(total) FILTER (WHERE value in (3, 4)) / sum(total), 2) AS pct_esp_vert_2015
    FROM cnt15
    WHERE value > 0 -- discard no data, including postgis raster no data
    GROUP BY "GeoUID"
),
cnt11 AS (
    SELECT "GeoUID"
        ,(pvc).value, SUM((pvc).count) As total
    FROM (SELECT "GeoUID"
            ,ST_ValueCount(ST_Clip(rast, geometry)) As pvc
        FROM canopee2011_10m
        JOIN ct16 ON ST_Intersects(geometry, rast)
    ) As foo
    GROUP BY "GeoUID", (pvc).value
),
canopee11 AS (
    SELECT "GeoUID"
        ,round(100. * sum(total) FILTER (WHERE value in (3, 4)) / sum(total), 2) AS pct_esp_vert_2011
    FROM cnt11
    WHERE value > 0 -- discard no data, including postgis raster no data
    GROUP BY "GeoUID"
)
SELECT "GeoUID"
    ,st_area(geometry) ct_area_m2
    ,coalesce(pct_esp_vert_2011, 0) pct_esp_vert_2011
    ,coalesce(pct_esp_vert_2015, 0) pct_esp_vert_2015
    ,coalesce(pct_esp_vert_2019, 0) pct_esp_vert_2019
FROM ct16
FULL JOIN canopee19 USING ("GeoUID")
FULL JOIN canopee15 USING ("GeoUID")
FULL JOIN canopee11 USING ("GeoUID");

1.2.1.3 Buffers 500m

WITH ct16 AS (
    select "GeoUID", ST_Buffer(geometry, 500) geometry
    from ct16
),
cnt19 AS (
    SELECT "GeoUID"
        ,(pvc).value, SUM((pvc).count) As total
    FROM (SELECT "GeoUID"
            ,ST_ValueCount(ST_Clip(rast, geometry)) As pvc
        FROM canopee2019_10m
        JOIN ct16 ON ST_Intersects(geometry, rast)
    ) As foo
    GROUP BY "GeoUID", (pvc).value
),
canopee19 AS (
    SELECT "GeoUID"
        ,round(100. * sum(total) FILTER (WHERE value in (3, 4)) / sum(total), 2) AS pct_esp_vert_2019
    FROM cnt19
    WHERE value > 0 -- discard no data, including postgis raster no data
    GROUP BY "GeoUID"
),
cnt15 AS (
    SELECT "GeoUID"
        ,(pvc).value, SUM((pvc).count) As total
    FROM (SELECT "GeoUID"
            ,ST_ValueCount(ST_Clip(rast, geometry)) As pvc
        FROM canopee2015_10m
        JOIN ct16 ON ST_Intersects(geometry, rast)
    ) As foo
    GROUP BY "GeoUID", (pvc).value
),
canopee15 AS (
    SELECT "GeoUID"
        ,round(100. * sum(total) FILTER (WHERE value in (3, 4)) / sum(total), 2) AS pct_esp_vert_2015
    FROM cnt15
    WHERE value > 0 -- discard no data, including postgis raster no data
    GROUP BY "GeoUID"
),
cnt11 AS (
    SELECT "GeoUID"
        ,(pvc).value, SUM((pvc).count) As total
    FROM (SELECT "GeoUID"
            ,ST_ValueCount(ST_Clip(rast, geometry)) As pvc
        FROM canopee2011_10m
        JOIN ct16 ON ST_Intersects(geometry, rast)
    ) As foo
    GROUP BY "GeoUID", (pvc).value
),
canopee11 AS (
    SELECT "GeoUID"
        ,round(100. * sum(total) FILTER (WHERE value in (3, 4)) / sum(total), 2) AS pct_esp_vert_2011
    FROM cnt11
    WHERE value > 0 -- discard no data, including postgis raster no data
    GROUP BY "GeoUID"
)
SELECT "GeoUID"
    ,st_area(geometry) ct_area_m2
    ,coalesce(pct_esp_vert_2011, 0) pct_esp_vert_2011
    ,coalesce(pct_esp_vert_2015, 0) pct_esp_vert_2015
    ,coalesce(pct_esp_vert_2019, 0) pct_esp_vert_2019
FROM ct16
FULL JOIN canopee19 USING ("GeoUID")
FULL JOIN canopee15 USING ("GeoUID")
FULL JOIN canopee11 USING ("GeoUID");

1.2.1.4 Buffers 750m

WITH ct16 AS (
    select "GeoUID", ST_Buffer(geometry, 750) geometry
    from ct16
),
cnt19 AS (
    SELECT "GeoUID"
        ,(pvc).value, SUM((pvc).count) As total
    FROM (SELECT "GeoUID"
            ,ST_ValueCount(ST_Clip(rast, geometry)) As pvc
        FROM canopee2019_10m
        JOIN ct16 ON ST_Intersects(geometry, rast)
    ) As foo
    GROUP BY "GeoUID", (pvc).value
),
canopee19 AS (
    SELECT "GeoUID"
        ,round(100. * sum(total) FILTER (WHERE value in (3, 4)) / sum(total), 2) AS pct_esp_vert_2019
    FROM cnt19
    WHERE value > 0 -- discard no data, including postgis raster no data
    GROUP BY "GeoUID"
),
cnt15 AS (
    SELECT "GeoUID"
        ,(pvc).value, SUM((pvc).count) As total
    FROM (SELECT "GeoUID"
            ,ST_ValueCount(ST_Clip(rast, geometry)) As pvc
        FROM canopee2015_10m
        JOIN ct16 ON ST_Intersects(geometry, rast)
    ) As foo
    GROUP BY "GeoUID", (pvc).value
),
canopee15 AS (
    SELECT "GeoUID"
        ,round(100. * sum(total) FILTER (WHERE value in (3, 4)) / sum(total), 2) AS pct_esp_vert_2015
    FROM cnt15
    WHERE value > 0 -- discard no data, including postgis raster no data
    GROUP BY "GeoUID"
),
cnt11 AS (
    SELECT "GeoUID"
        ,(pvc).value, SUM((pvc).count) As total
    FROM (SELECT "GeoUID"
            ,ST_ValueCount(ST_Clip(rast, geometry)) As pvc
        FROM canopee2011_10m
        JOIN ct16 ON ST_Intersects(geometry, rast)
    ) As foo
    GROUP BY "GeoUID", (pvc).value
),
canopee11 AS (
    SELECT "GeoUID"
        ,round(100. * sum(total) FILTER (WHERE value in (3, 4)) / sum(total), 2) AS pct_esp_vert_2011
    FROM cnt11
    WHERE value > 0 -- discard no data, including postgis raster no data
    GROUP BY "GeoUID"
)
SELECT "GeoUID"
    ,st_area(geometry) ct_area_m2
    ,coalesce(pct_esp_vert_2011, 0) pct_esp_vert_2011
    ,coalesce(pct_esp_vert_2015, 0) pct_esp_vert_2015
    ,coalesce(pct_esp_vert_2019, 0) pct_esp_vert_2019
FROM ct16
FULL JOIN canopee19 USING ("GeoUID")
FULL JOIN canopee15 USING ("GeoUID")
FULL JOIN canopee11 USING ("GeoUID");

1.3 Pampalon index

Get it here

pampalon <- read.xlsx("data/Canada2016Pampalon/A-MSDIData_Can2016_eng/1. EquivalenceTableCanada2016_ENG.xlsx", sheet = 2) %>%
  mutate(DA = as.character(DA)) %>%
  select(DA, SCOREMAT, SCORESOC)

# 2016 DA boundaries for Montreal
DA16 <- get_census(dataset='CA16', regions=list(CMA='24462'), level='DA', geo_format = "sf") %>%
  filter(Type == "DA") %>%
  st_transform(st_crs(bike_lane))

## Reading geo data from local cache.

pampalon <- DA16 %>%
  inner_join(pampalon, by = c("GeoUID" = "DA")) %>%
  as.data.frame()

# Get Pampalon 2011
pampalon11 <- read.xlsx("data/Canada2011Pampalon/A-MSDIData_Can2011_eng/1. CorrespondenceTable_Can2011_eng.xlsx", sheet = 2) %>%
  mutate(DA = as.character(DA)) %>%
  select(DA, DAPOP2011, SCOREMAT, SCORESOC)

# Get LUT DA2011 <-> DA2016 from StatCan
lut_da <- read.csv("data/2016_92-156_DA_AD_csv/2016_92-156_DA_AD.csv", colClasses = "character") %>%
  select(!c(DBUID2016.IDIDU2016, DA_rel_flag.AD_ind_rel)) %>%
  unique()

# Link Pampalon 2011 to LUT and compute weighted mean of scores of Pampalon 2011
pampalon11.16 <- pampalon11 %>%
  inner_join(lut_da, by = c("DA" = "DAUID2011.ADIDU2011")) %>%
  group_by(DAUID2016.ADIDU2016) %>%
  summarise(pop2011 = sum(DAPOP2011),
            SCOREMAT = weighted.mean(SCOREMAT, DAPOP2011, na.rm = TRUE),
            SCORESOC = weighted.mean(SCORESOC, DAPOP2011, na.rm = TRUE))

# Then link Pampalon 2011 to 2016
pampalon <- pampalon %>%
  left_join(pampalon11.16, by = c("GeoUID" = "DAUID2016.ADIDU2016"), suffix = c(".16", ".11"))

# Aggregate at the CT level
pampalon_CT <- pampalon %>%
  group_by(CT_UID) %>%
  summarise(wSCOREMAT.2016 = weighted.mean(SCOREMAT.16, Population, na.rm = TRUE),
            wSCORESOC.2016 = weighted.mean(SCORESOC.16, Population, na.rm = TRUE),
            wSCOREMAT.2011 = weighted.mean(SCOREMAT.11, pop2011, na.rm = TRUE),
            wSCORESOC.2011 = weighted.mean(SCORESOC.11, pop2011, na.rm = TRUE))

# Clean up
rm(lut_da, pampalon11.16, pampalon11)

2 Preliminary analyses

Outcome variables = Pampalon deprivation index or gentrification flag Independent variables = all bike lane and canopy variables

2.1 Build complete dataset

All variables + outcome linked at the CT level

.bike_lane_changes <- bike_lane_changes %>%
  as.data.frame() %>%
  select(GeoUID, ends_with("ct", ignore.case = FALSE), ends_with("b250", ignore.case = FALSE), ends_with("b500", ignore.case = FALSE), ends_with("b750", ignore.case = FALSE)) %>%
  select(GeoUID, starts_with("Bike_lane")) # Drop individual category lane length

bei_df <- CT16 %>%
  as.data.frame() %>%
  transmute(CT_UID = GeoUID,
            CD_UID = CD_UID,
            CSD_UID = CSD_UID,
            interact_aoi = interact_aoi,
            Population = Population) %>%
  left_join(pampalon_CT, by="CT_UID") %>%
  left_join(.bike_lane_changes, by=c("CT_UID" = "GeoUID")) %>%
  left_join(as.data.frame(esp_vert_ct), by=c("CT_UID" = "GeoUID")) %>%
  left_join(as.data.frame(esp_vert_buf250), by=c("CT_UID" = "GeoUID"), suffix=c("ct", "b250")) %>%
  left_join(as.data.frame(esp_vert_buf500), by=c("CT_UID" = "GeoUID")) %>%
  left_join(as.data.frame(esp_vert_buf750), by=c("CT_UID" = "GeoUID"), suffix=c("b500", "b750")) %>%
  left_join(select(as.data.frame(ding$`2016`), ct_uid_16, starts_with("gentrif")), by=c("CT_UID" = "ct_uid_16")) %>%
  left_join(select(as.data.frame(ding$`2011`), ct_uid_11, starts_with("gentrif")), by=c("CT_UID" = "ct_uid_11")) %>%
  left_join(select(as.data.frame(ding$`2006`), ct_uid_06, starts_with("gentrif")), by=c("CT_UID" = "ct_uid_06"))
  
  
head(bei_df)

write.csv(bei_df, "data/_results/bei_equity.csv", na="", row.names = FALSE)

Included variables:

• Census Tracts variables
  • CT_UID: 2016 Census Tract ID
  • CD_UID: 2016 Census Division
  • CSD_UID: 2016 Census Subdivision
  • interact_aoi: Does CT belong to INTERACT study area?
  • Population: 2016 Population within CT
  • ct_area_m2.{ct|b{250|500|750}}: Area of CT or buffer of 250, 500 or 750m radius around CT, in square meters
  • gentrified_2016_2011: Is the CT gentrified in 2016?
  • gentrifiable_2011: Is the CT candidate to gentrification in 2011?
  • gentrified_2011_2006: Is the CT gentrified in 2011
  • gentrifiable_2006: Is the CT candidate to gentrification in 2006
  • gentrified_2006_2001: Is the CT gentrified in 2006
  • gentrifiable_2001: Is the CT candidate to gentrification in 2001
• Pampalon’s metrics
  • wSCOREMAT.2016: Social deprivation index in 2016 (population weighted)
  • wSCORESOC.2016: Material deprivation index in 2016 (population weighted)
  • wSCOREMAT.2011: Social deprivation index in 2011 (population weighted)
  • wSCORESOC.2011: Material deprivation index in 2011 (population weighted)
• Bike lane density
  • Bike_lane_total.{2016|2011}{ct|b{250|500|750}}: total length of bike lanes, in 2016 or 2011, within CT or buffer of 250, 500 or 750m radius
  • Bike_lane_diff.2011.2016{ct|b{250|500|750}}: change in total length of bike lanes between 2011 and 2011, within CT or buffer of 250, 500 or 750m radius
  • Bike_lane_diff.by.street.2011.2016{ct|b{250|500|750}}: change in total length of bike lanes between 2011 and 2011, normalized by street length, within CT or buffer of 250, 500 or 750m radius
  • Bike_lane_diff.by.area.2011.2016{ct|b{250|500|750}}: change in total length of bike lanes between 2011 and 2011, normalized by area, within CT or buffer of 250, 500 or 750m radius
• Green spaces
  • pct_esp_vert_{2011|2015|2019}.{ct|b{250|500|750}}: % of green space in 2011, 2015 or 2019 within CT or buffer of 250, 500 or 750m radius
  • pct_esp_vert_diff{2011|2015}.{2015|2019}.{ct|b{250|500|750}}: change in % of green space between 2011 and 2015, 2011 and 2019 as well as 2011 and 2019, within CT or buffer of 250, 500 or 750m radius

2.2 Preliminary analyses

2.2.1 Outcome variable distribution

Whole area ~ Montréal CMA / CMM

.bei_df_long <- bei_df %>% 
  units::drop_units() %>% 
  select(CT_UID, CD_UID, Population, starts_with("wSCORE"), starts_with("gentrif")) %>%
  pivot_longer(!c(CT_UID, CD_UID))

ggplot(.bei_df_long, aes(value)) +
  geom_histogram() + 
  facet_wrap(~name, scales = "free", ncol = 3)

## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

## Warning: Removed 618 rows containing non-finite values (stat_bin).

2.2.2 Exposure variable distributions

2.2.2.1 Census Tracts

.bei_df_long <- bei_df %>% 
  units::drop_units() %>% 
  select(CT_UID, CD_UID, matches("^Bike_lane.*ct$"), matches("^pct_esp_vert_diff.*ct$")) %>%
  pivot_longer(!c(CT_UID, CD_UID))

ggplot(.bei_df_long, aes(value)) +
  geom_histogram() + 
  facet_wrap(~name, scales = "free")

## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

## Warning: Removed 252 rows containing non-finite values (stat_bin).

2.2.2.2 Buffers 250m

.bei_df_long <- bei_df %>% 
  units::drop_units() %>% 
  select(CT_UID, CD_UID, matches("^Bike_lane.*b250$"), matches("^pct_esp_vert_diff.*b250$")) %>%
  pivot_longer(!c(CT_UID, CD_UID))

ggplot(.bei_df_long, aes(value)) +
  geom_histogram() + 
  facet_wrap(~name, scales = "free")

## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

## Warning: Removed 233 rows containing non-finite values (stat_bin).

2.2.2.3 Buffers 500m

.bei_df_long <- bei_df %>% 
  units::drop_units() %>% 
  select(CT_UID, CD_UID, matches("^Bike_lane.*b500$"), matches("^pct_esp_vert_diff.*b500$")) %>%
  pivot_longer(!c(CT_UID, CD_UID))

ggplot(.bei_df_long, aes(value)) +
  geom_histogram() + 
  facet_wrap(~name, scales = "free")

## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

## Warning: Removed 229 rows containing non-finite values (stat_bin).

2.2.2.4 Buffers 750m

.bei_df_long <- bei_df %>% 
  units::drop_units() %>% 
  select(CT_UID, CD_UID, matches("^Bike_lane.*b750$"), matches("^pct_esp_vert_diff.*b750$")) %>%
  pivot_longer(!c(CT_UID, CD_UID))

ggplot(.bei_df_long, aes(value)) +
  geom_histogram() + 
  facet_wrap(~name, scales = "free")

## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

## Warning: Removed 222 rows containing non-finite values (stat_bin).

2.2.3 Outcome variable distribution | INTERACT study area

INTERACT study area ~ Montréal, Laval, Longueuil, Brossard, St-Lambert

.bei_df_long <- bei_df %>% 
  filter(interact_aoi) %>%
  units::drop_units() %>% 
  select(CT_UID, CD_UID, Population, starts_with("wSCORE"), starts_with("gentrif")) %>%
  pivot_longer(!c(CT_UID, CD_UID))

ggplot(.bei_df_long, aes(value)) +
  geom_histogram() + 
  facet_wrap(~name, scales = "free", ncol = 3)

## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

## Warning: Removed 276 rows containing non-finite values (stat_bin).

2.2.4 Exposure variable distributions | INTERACT study area

2.2.4.1 Census Tracts

.bei_df_long <- bei_df %>% 
  filter(interact_aoi) %>%
  units::drop_units() %>% 
  select(CT_UID, CD_UID, matches("^Bike_lane.*ct$"), matches("^pct_esp_vert_diff.*ct$")) %>%
  pivot_longer(!c(CT_UID, CD_UID))

ggplot(.bei_df_long, aes(value)) +
  geom_histogram() + 
  facet_wrap(~name, scales = "free")

## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

2.2.4.2 Buffers 250m

.bei_df_long <- bei_df %>% 
  filter(interact_aoi) %>%
  units::drop_units() %>% 
  select(CT_UID, CD_UID, matches("^Bike_lane.*b250$"), matches("^pct_esp_vert_diff.*b250$")) %>%
  pivot_longer(!c(CT_UID, CD_UID))

ggplot(.bei_df_long, aes(value)) +
  geom_histogram() + 
  facet_wrap(~name, scales = "free")

## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

2.2.4.3 Buffers 500m

.bei_df_long <- bei_df %>% 
  filter(interact_aoi) %>%
  units::drop_units() %>% 
  select(CT_UID, CD_UID, matches("^Bike_lane.*b500$"), matches("^pct_esp_vert_diff.*b500$")) %>%
  pivot_longer(!c(CT_UID, CD_UID))

ggplot(.bei_df_long, aes(value)) +
  geom_histogram() + 
  facet_wrap(~name, scales = "free")

## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

2.2.4.4 Buffers 750m

.bei_df_long <- bei_df %>% 
  filter(interact_aoi) %>%
  units::drop_units() %>% 
  select(CT_UID, CD_UID, matches("^Bike_lane.*b750$"), matches("^pct_esp_vert_diff.*b750$")) %>%
  pivot_longer(!c(CT_UID, CD_UID))

ggplot(.bei_df_long, aes(value)) +
  geom_histogram() + 
  facet_wrap(~name, scales = "free")

## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

2.3 Association between outcomes and independant variables

Questions @Yan:

• veut-on synchroniser les périodes de gentrification avec l’exposition sur le terrain (pour l’instant on se contente de regarder les différentes périodes de gentrification avec la même période d’observation des changements d’exposition – bike lanes changes 2011 -> 2016 / greenspace changes 2011 -> 2019)
• veut-on changer la fenêtre de mesure des changements ? (pour l’instant, 5 ans pour gentrif’ & bike lane changes / 8 ans pour greenness)

2.3.1 Pampalon outcome

2.3.1.1 Census Tracts

2.3.1.1.1 Pampalon / Material vs. bike lane length change

ggplot(units::drop_units(bei_df), aes(x=Bike_lane_diff.2011.2016ct, y=wSCOREMAT.2016)) +
  geom_point() +
  geom_smooth(method=lm)

res.lm <- lm(wSCOREMAT.2016 ~ Bike_lane_diff.2011.2016ct, data = units::drop_units(bei_df))
summary(res.lm)

## 
## Call:
## lm(formula = wSCOREMAT.2016 ~ Bike_lane_diff.2011.2016ct, data = units::drop_units(bei_df))
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.14093 -0.02165 -0.00058  0.02353  0.13849 
## 
## Coefficients:
##                              Estimate Std. Error t value Pr(>|t|)  
## (Intercept)                 3.904e-04  1.219e-03    0.32   0.7488  
## Bike_lane_diff.2011.2016ct -3.201e-06  1.301e-06   -2.46   0.0141 *
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.03602 on 948 degrees of freedom
##   (20 observations deleted due to missingness)
## Multiple R-squared:  0.006345,   Adjusted R-squared:  0.005297 
## F-statistic: 6.053 on 1 and 948 DF,  p-value: 0.01406

2.3.1.1.2 Pampalon / Social vs. bike lane length change

ggplot(units::drop_units(bei_df), aes(x=Bike_lane_diff.2011.2016ct, y=wSCORESOC.2016)) +
  geom_point() +
  geom_smooth(method=lm)

res.lm <- lm(wSCORESOC.2016 ~ Bike_lane_diff.2011.2016ct, data = units::drop_units(bei_df))
summary(res.lm)

## 
## Call:
## lm(formula = wSCORESOC.2016 ~ Bike_lane_diff.2011.2016ct, data = units::drop_units(bei_df))
## 
## Residuals:
##       Min        1Q    Median        3Q       Max 
## -0.076895 -0.022870  0.000003  0.024011  0.094864 
## 
## Coefficients:
##                              Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                 1.810e-02  9.833e-04  18.411   <2e-16 ***
## Bike_lane_diff.2011.2016ct -2.671e-07  1.050e-06  -0.254    0.799    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.02906 on 948 degrees of freedom
##   (20 observations deleted due to missingness)
## Multiple R-squared:  6.83e-05,   Adjusted R-squared:  -0.0009865 
## F-statistic: 0.06475 on 1 and 948 DF,  p-value: 0.7992

2.3.1.1.3 Pampalon / Material vs. bike lane length change, normalized by street length

ggplot(units::drop_units(bei_df), aes(x=Bike_lane_diff.by.street.2011.2016ct, y=wSCOREMAT.2016)) +
  geom_point() +
  geom_smooth(method=lm)

res.lm <- lm(wSCOREMAT.2016 ~ Bike_lane_diff.by.street.2011.2016ct, data = units::drop_units(bei_df))
summary(res.lm)

## 
## Call:
## lm(formula = wSCOREMAT.2016 ~ Bike_lane_diff.by.street.2011.2016ct, 
##     data = units::drop_units(bei_df))
## 
## Residuals:
##       Min        1Q    Median        3Q       Max 
## -0.140946 -0.025988  0.001677  0.026813  0.138472 
## 
## Coefficients:
##                                        Estimate Std. Error t value Pr(>|t|)
## (Intercept)                           0.0004043  0.0016064   0.252    0.801
## Bike_lane_diff.by.street.2011.2016ct -0.0164586  0.0155086  -1.061    0.289
## 
## Residual standard error: 0.03978 on 702 degrees of freedom
##   (266 observations deleted due to missingness)
## Multiple R-squared:  0.001602,   Adjusted R-squared:  0.0001796 
## F-statistic: 1.126 on 1 and 702 DF,  p-value: 0.2889

2.3.1.1.4 Pampalon / Social vs. bike lane length change, normalized by street length

ggplot(units::drop_units(bei_df), aes(x=Bike_lane_diff.by.street.2011.2016ct, y=wSCORESOC.2016)) +
  geom_point() +
  geom_smooth(method=lm)

res.lm <- lm(wSCORESOC.2016 ~ Bike_lane_diff.by.street.2011.2016ct, data = units::drop_units(bei_df))
summary(res.lm)

## 
## Call:
## lm(formula = wSCORESOC.2016 ~ Bike_lane_diff.by.street.2011.2016ct, 
##     data = units::drop_units(bei_df))
## 
## Residuals:
##       Min        1Q    Median        3Q       Max 
## -0.099465 -0.020801  0.003917  0.021838  0.091426 
## 
## Coefficients:
##                                      Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                          0.021544   0.001169  18.423  < 2e-16 ***
## Bike_lane_diff.by.street.2011.2016ct 0.038459   0.011290   3.407 0.000695 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.02896 on 702 degrees of freedom
##   (266 observations deleted due to missingness)
## Multiple R-squared:  0.01626,    Adjusted R-squared:  0.01486 
## F-statistic:  11.6 on 1 and 702 DF,  p-value: 0.0006954

2.3.1.1.5 Pampalon / Material vs. bike lane length change, normalized by area

ggplot(units::drop_units(bei_df), aes(x=Bike_lane_diff.by.area.2011.2016ct, y=wSCOREMAT.2016)) +
  geom_point() +
  geom_smooth(method=lm)

res.lm <- lm(wSCOREMAT.2016 ~ Bike_lane_diff.by.area.2011.2016ct, data = units::drop_units(bei_df))
summary(res.lm)

## 
## Call:
## lm(formula = wSCOREMAT.2016 ~ Bike_lane_diff.by.area.2011.2016ct, 
##     data = units::drop_units(bei_df))
## 
## Residuals:
##       Min        1Q    Median        3Q       Max 
## -0.140461 -0.022162 -0.000097  0.023480  0.138958 
## 
## Coefficients:
##                                      Estimate Std. Error t value Pr(>|t|)
## (Intercept)                        -8.102e-05  1.265e-03  -0.064    0.949
## Bike_lane_diff.by.area.2011.2016ct -8.785e-04  1.104e-03  -0.796    0.426
## 
## Residual standard error: 0.03612 on 948 degrees of freedom
##   (20 observations deleted due to missingness)
## Multiple R-squared:  0.0006678,  Adjusted R-squared:  -0.0003864 
## F-statistic: 0.6335 on 1 and 948 DF,  p-value: 0.4263

2.3.1.1.6 Pampalon / Social vs. bike lane length change, normalized by area

ggplot(units::drop_units(bei_df), aes(x=Bike_lane_diff.by.area.2011.2016ct, y=wSCORESOC.2016)) +
  geom_point() +
  geom_smooth(method=lm)

res.lm <- lm(wSCORESOC.2016 ~ Bike_lane_diff.by.area.2011.2016ct, data = units::drop_units(bei_df))
summary(res.lm)

## 
## Call:
## lm(formula = wSCORESOC.2016 ~ Bike_lane_diff.by.area.2011.2016ct, 
##     data = units::drop_units(bei_df))
## 
## Residuals:
##       Min        1Q    Median        3Q       Max 
## -0.073460 -0.021104  0.000026  0.019814  0.098300 
## 
## Coefficients:
##                                     Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                        0.0146689  0.0009759  15.031   <2e-16 ***
## Bike_lane_diff.by.area.2011.2016ct 0.0077844  0.0008513   9.144   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.02786 on 948 degrees of freedom
##   (20 observations deleted due to missingness)
## Multiple R-squared:  0.08105,    Adjusted R-squared:  0.08008 
## F-statistic: 83.61 on 1 and 948 DF,  p-value: < 2.2e-16

2.3.1.1.7 Pampalon / Material vs. green space change

ggplot(units::drop_units(bei_df), aes(x=pct_esp_vert_diff_2011.2019ct, y=wSCOREMAT.2016)) +
  geom_point() +
  geom_smooth(method=lm)

res.lm <- lm(wSCOREMAT.2016 ~ pct_esp_vert_diff_2011.2019ct, data = units::drop_units(bei_df))
summary(res.lm)

## 
## Call:
## lm(formula = wSCOREMAT.2016 ~ pct_esp_vert_diff_2011.2019ct, 
##     data = units::drop_units(bei_df))
## 
## Residuals:
##       Min        1Q    Median        3Q       Max 
## -0.140498 -0.022107  0.000299  0.023728  0.142023 
## 
## Coefficients:
##                                 Estimate Std. Error t value Pr(>|t|)
## (Intercept)                    0.0004029  0.0014079   0.286    0.775
## pct_esp_vert_diff_2011.2019ct -0.0002172  0.0001964  -1.106    0.269
## 
## Residual standard error: 0.03611 on 948 degrees of freedom
##   (20 observations deleted due to missingness)
## Multiple R-squared:  0.001289,   Adjusted R-squared:  0.0002352 
## F-statistic: 1.223 on 1 and 948 DF,  p-value: 0.269

2.3.1.1.8 Pampalon / Social vs. green space change

ggplot(units::drop_units(bei_df), aes(x=pct_esp_vert_diff_2011.2019ct, y=wSCORESOC.2016)) +
  geom_point() +
  geom_smooth(method=lm)

res.lm <- lm(wSCORESOC.2016 ~ pct_esp_vert_diff_2011.2019ct, data = units::drop_units(bei_df))
summary(res.lm)

## 
## Call:
## lm(formula = wSCORESOC.2016 ~ pct_esp_vert_diff_2011.2019ct, 
##     data = units::drop_units(bei_df))
## 
## Residuals:
##       Min        1Q    Median        3Q       Max 
## -0.076679 -0.022820 -0.000237  0.024013  0.095232 
## 
## Coefficients:
##                                Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                   0.0174185  0.0011326  15.379   <2e-16 ***
## pct_esp_vert_diff_2011.2019ct 0.0001547  0.0001580   0.979    0.328    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.02905 on 948 degrees of freedom
##   (20 observations deleted due to missingness)
## Multiple R-squared:  0.00101,    Adjusted R-squared:  -4.349e-05 
## F-statistic: 0.9587 on 1 and 948 DF,  p-value: 0.3278

2.3.1.2 Buffers 250m

2.3.1.2.1 Pampalon / Material vs. bike lane length change

ggplot(units::drop_units(bei_df), aes(x=Bike_lane_diff.2011.2016b250, y=wSCOREMAT.2016)) +
  geom_point() +
  geom_smooth(method=lm)

res.lm <- lm(wSCOREMAT.2016 ~ Bike_lane_diff.2011.2016b250, data = units::drop_units(bei_df))
summary(res.lm)

## 
## Call:
## lm(formula = wSCOREMAT.2016 ~ Bike_lane_diff.2011.2016b250, data = units::drop_units(bei_df))
## 
## Residuals:
##       Min        1Q    Median        3Q       Max 
## -0.141558 -0.021775 -0.000925  0.023297  0.137860 
## 
## Coefficients:
##                                Estimate Std. Error t value Pr(>|t|)   
## (Intercept)                   1.016e-03  1.277e-03   0.796  0.42622   
## Bike_lane_diff.2011.2016b250 -2.258e-06  7.906e-07  -2.856  0.00438 **
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.03598 on 948 degrees of freedom
##   (20 observations deleted due to missingness)
## Multiple R-squared:  0.008533,   Adjusted R-squared:  0.007487 
## F-statistic: 8.159 on 1 and 948 DF,  p-value: 0.004378

2.3.1.2.2 Pampalon / Social vs. bike lane length change

ggplot(units::drop_units(bei_df), aes(x=Bike_lane_diff.2011.2016b250, y=wSCORESOC.2016)) +
  geom_point() +
  geom_smooth(method=lm)

res.lm <- lm(wSCORESOC.2016 ~ Bike_lane_diff.2011.2016b250, data = units::drop_units(bei_df))
summary(res.lm)

## 
## Call:
## lm(formula = wSCORESOC.2016 ~ Bike_lane_diff.2011.2016b250, data = units::drop_units(bei_df))
## 
## Residuals:
##       Min        1Q    Median        3Q       Max 
## -0.075420 -0.022369  0.000121  0.023588  0.096340 
## 
## Coefficients:
##                               Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                  1.663e-02  1.025e-03  16.221  < 2e-16 ***
## Bike_lane_diff.2011.2016b250 2.148e-06  6.348e-07   3.383 0.000746 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.02889 on 948 degrees of freedom
##   (20 observations deleted due to missingness)
## Multiple R-squared:  0.01193,    Adjusted R-squared:  0.01089 
## F-statistic: 11.45 on 1 and 948 DF,  p-value: 0.0007459

2.3.1.2.3 Pampalon / Material vs. bike lane length change, normalized by street length

ggplot(units::drop_units(bei_df), aes(x=Bike_lane_diff.by.street.2011.2016b250, y=wSCOREMAT.2016)) +
  geom_point() +
  geom_smooth(method=lm)

res.lm <- lm(wSCOREMAT.2016 ~ Bike_lane_diff.by.street.2011.2016b250, data = units::drop_units(bei_df))
summary(res.lm)

## 
## Call:
## lm(formula = wSCOREMAT.2016 ~ Bike_lane_diff.by.street.2011.2016b250, 
##     data = units::drop_units(bei_df))
## 
## Residuals:
##       Min        1Q    Median        3Q       Max 
## -0.142037 -0.025277  0.001961  0.026047  0.137382 
## 
## Coefficients:
##                                         Estimate Std. Error t value Pr(>|t|)  
## (Intercept)                             0.001495   0.001676   0.892    0.373  
## Bike_lane_diff.by.street.2011.2016b250 -0.055109   0.024530  -2.247    0.025 *
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.03926 on 721 degrees of freedom
##   (247 observations deleted due to missingness)
## Multiple R-squared:  0.006952,   Adjusted R-squared:  0.005574 
## F-statistic: 5.047 on 1 and 721 DF,  p-value: 0.02497

2.3.1.2.4 Pampalon / Social vs. bike lane length change, normalized by street length

ggplot(units::drop_units(bei_df), aes(x=Bike_lane_diff.by.street.2011.2016b250, y=wSCORESOC.2016)) +
  geom_point() +
  geom_smooth(method=lm)

res.lm <- lm(wSCORESOC.2016 ~ Bike_lane_diff.by.street.2011.2016b250, data = units::drop_units(bei_df))
summary(res.lm)

## 
## Call:
## lm(formula = wSCORESOC.2016 ~ Bike_lane_diff.by.street.2011.2016b250, 
##     data = units::drop_units(bei_df))
## 
## Residuals:
##       Min        1Q    Median        3Q       Max 
## -0.085760 -0.020329  0.003543  0.020790  0.094580 
## 
## Coefficients:
##                                        Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                            0.018389   0.001209  15.210  < 2e-16 ***
## Bike_lane_diff.by.street.2011.2016b250 0.124323   0.017693   7.027 4.91e-12 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.02832 on 721 degrees of freedom
##   (247 observations deleted due to missingness)
## Multiple R-squared:  0.06409,    Adjusted R-squared:  0.06279 
## F-statistic: 49.37 on 1 and 721 DF,  p-value: 4.906e-12

2.3.1.2.5 Pampalon / Material vs. bike lane length change, normalized by area

ggplot(units::drop_units(bei_df), aes(x=Bike_lane_diff.by.area.2011.2016b250, y=wSCOREMAT.2016)) +
  geom_point() +
  geom_smooth(method=lm)

res.lm <- lm(wSCOREMAT.2016 ~ Bike_lane_diff.by.area.2011.2016b250, data = units::drop_units(bei_df))
summary(res.lm)

## 
## Call:
## lm(formula = wSCOREMAT.2016 ~ Bike_lane_diff.by.area.2011.2016b250, 
##     data = units::drop_units(bei_df))
## 
## Residuals:
##       Min        1Q    Median        3Q       Max 
## -0.141396 -0.022428 -0.000406  0.023373  0.138023 
## 
## Coefficients:
##                                        Estimate Std. Error t value Pr(>|t|)  
## (Intercept)                           0.0008536  0.0013002   0.657   0.5117  
## Bike_lane_diff.by.area.2011.2016b250 -0.0032359  0.0014014  -2.309   0.0212 *
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.03603 on 948 degrees of freedom
##   (20 observations deleted due to missingness)
## Multiple R-squared:  0.005593,   Adjusted R-squared:  0.004544 
## F-statistic: 5.332 on 1 and 948 DF,  p-value: 0.02115

2.3.1.2.6 Pampalon / Social vs. bike lane length change, normalized by area

ggplot(units::drop_units(bei_df), aes(x=Bike_lane_diff.by.area.2011.2016b250, y=wSCORESOC.2016)) +
  geom_point() +
  geom_smooth(method=lm)

res.lm <- lm(wSCORESOC.2016 ~ Bike_lane_diff.by.area.2011.2016b250, data = units::drop_units(bei_df))
summary(res.lm)

## 
## Call:
## lm(formula = wSCORESOC.2016 ~ Bike_lane_diff.by.area.2011.2016b250, 
##     data = units::drop_units(bei_df))
## 
## Residuals:
##       Min        1Q    Median        3Q       Max 
## -0.071772 -0.019890 -0.000031  0.019665  0.099988 
## 
## Coefficients:
##                                       Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                          0.0129809  0.0009795   13.25   <2e-16 ***
## Bike_lane_diff.by.area.2011.2016b250 0.0124400  0.0010557   11.78   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.02714 on 948 degrees of freedom
##   (20 observations deleted due to missingness)
## Multiple R-squared:  0.1278, Adjusted R-squared:  0.1268 
## F-statistic: 138.8 on 1 and 948 DF,  p-value: < 2.2e-16

2.3.1.2.7 Pampalon / Material vs. green space change

ggplot(units::drop_units(bei_df), aes(x=pct_esp_vert_diff_2011.2019b250, y=wSCOREMAT.2016)) +
  geom_point() +
  geom_smooth(method=lm)

res.lm <- lm(wSCOREMAT.2016 ~ pct_esp_vert_diff_2011.2019b250, data = units::drop_units(bei_df))
summary(res.lm)

## 
## Call:
## lm(formula = wSCOREMAT.2016 ~ pct_esp_vert_diff_2011.2019b250, 
##     data = units::drop_units(bei_df))
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.14025 -0.02209  0.00034  0.02352  0.14000 
## 
## Coefficients:
##                                   Estimate Std. Error t value Pr(>|t|)  
## (Intercept)                      0.0009586  0.0014049   0.682   0.4952  
## pct_esp_vert_diff_2011.2019b250 -0.0003734  0.0002045  -1.826   0.0682 .
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.03607 on 948 degrees of freedom
##   (20 observations deleted due to missingness)
## Multiple R-squared:  0.003504,   Adjusted R-squared:  0.002453 
## F-statistic: 3.333 on 1 and 948 DF,  p-value: 0.0682

2.3.1.2.8 Pampalon / Social vs. green space change

ggplot(units::drop_units(bei_df), aes(x=pct_esp_vert_diff_2011.2019b250, y=wSCORESOC.2016)) +
  geom_point() +
  geom_smooth(method=lm)

res.lm <- lm(wSCORESOC.2016 ~ pct_esp_vert_diff_2011.2019b250, data = units::drop_units(bei_df))
summary(res.lm)

## 
## Call:
## lm(formula = wSCORESOC.2016 ~ pct_esp_vert_diff_2011.2019b250, 
##     data = units::drop_units(bei_df))
## 
## Residuals:
##       Min        1Q    Median        3Q       Max 
## -0.077329 -0.022934 -0.000117  0.024017  0.094963 
## 
## Coefficients:
##                                  Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                     1.781e-02  1.132e-03  15.731   <2e-16 ***
## pct_esp_vert_diff_2011.2019b250 5.969e-05  1.648e-04   0.362    0.717    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.02906 on 948 degrees of freedom
##   (20 observations deleted due to missingness)
## Multiple R-squared:  0.0001384,  Adjusted R-squared:  -0.0009163 
## F-statistic: 0.1312 on 1 and 948 DF,  p-value: 0.7173

2.3.1.3 Buffers 500m

2.3.1.3.1 Pampalon / Material vs. bike lane length change

ggplot(units::drop_units(bei_df), aes(x=Bike_lane_diff.2011.2016b500, y=wSCOREMAT.2016)) +
  geom_point() +
  geom_smooth(method=lm)

res.lm <- lm(wSCOREMAT.2016 ~ Bike_lane_diff.2011.2016b500, data = units::drop_units(bei_df))
summary(res.lm)

## 
## Call:
## lm(formula = wSCOREMAT.2016 ~ Bike_lane_diff.2011.2016b500, data = units::drop_units(bei_df))
## 
## Residuals:
##       Min        1Q    Median        3Q       Max 
## -0.142173 -0.021775 -0.000562  0.023345  0.137246 
## 
## Coefficients:
##                                Estimate Std. Error t value Pr(>|t|)   
## (Intercept)                   1.630e-03  1.328e-03   1.228  0.21993   
## Bike_lane_diff.2011.2016b500 -1.772e-06  5.396e-07  -3.284  0.00106 **
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.03593 on 948 degrees of freedom
##   (20 observations deleted due to missingness)
## Multiple R-squared:  0.01125,    Adjusted R-squared:  0.01021 
## F-statistic: 10.79 on 1 and 948 DF,  p-value: 0.00106

2.3.1.3.2 Pampalon / Social vs. bike lane length change

ggplot(units::drop_units(bei_df), aes(x=Bike_lane_diff.2011.2016b500, y=wSCORESOC.2016)) +
  geom_point() +
  geom_smooth(method=lm)

res.lm <- lm(wSCORESOC.2016 ~ Bike_lane_diff.2011.2016b500, data = units::drop_units(bei_df))
summary(res.lm)

## 
## Call:
## lm(formula = wSCORESOC.2016 ~ Bike_lane_diff.2011.2016b500, data = units::drop_units(bei_df))
## 
## Residuals:
##       Min        1Q    Median        3Q       Max 
## -0.087146 -0.021312  0.000382  0.021537  0.098092 
## 
## Coefficients:
##                               Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                  1.488e-02  1.053e-03  14.130  < 2e-16 ***
## Bike_lane_diff.2011.2016b500 2.675e-06  4.278e-07   6.253 6.09e-10 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.02848 on 948 degrees of freedom
##   (20 observations deleted due to missingness)
## Multiple R-squared:  0.03961,    Adjusted R-squared:  0.0386 
## F-statistic:  39.1 on 1 and 948 DF,  p-value: 6.089e-10

2.3.1.3.3 Pampalon / Material vs. green space change

ggplot(units::drop_units(bei_df), aes(x=pct_esp_vert_diff_2011.2019b500, y=wSCOREMAT.2016)) +
  geom_point() +
  geom_smooth(method=lm)

res.lm <- lm(wSCOREMAT.2016 ~ pct_esp_vert_diff_2011.2019b500, data = units::drop_units(bei_df))
summary(res.lm)

## 
## Call:
## lm(formula = wSCOREMAT.2016 ~ pct_esp_vert_diff_2011.2019b500, 
##     data = units::drop_units(bei_df))
## 
## Residuals:
##       Min        1Q    Median        3Q       Max 
## -0.139854 -0.022182  0.000328  0.023299  0.140003 
## 
## Coefficients:
##                                   Estimate Std. Error t value Pr(>|t|)  
## (Intercept)                      0.0009664  0.0014117   0.685    0.494  
## pct_esp_vert_diff_2011.2019b500 -0.0003804  0.0002105  -1.807    0.071 .
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.03607 on 948 degrees of freedom
##   (20 observations deleted due to missingness)
## Multiple R-squared:  0.003434,   Adjusted R-squared:  0.002383 
## F-statistic: 3.267 on 1 and 948 DF,  p-value: 0.07101

2.3.1.3.4 Pampalon / Social vs. green space change

ggplot(units::drop_units(bei_df), aes(x=pct_esp_vert_diff_2011.2019b500, y=wSCORESOC.2016)) +
  geom_point() +
  geom_smooth(method=lm)

res.lm <- lm(wSCORESOC.2016 ~ pct_esp_vert_diff_2011.2019b500, data = units::drop_units(bei_df))
summary(res.lm)

## 
## Call:
## lm(formula = wSCORESOC.2016 ~ pct_esp_vert_diff_2011.2019b500, 
##     data = units::drop_units(bei_df))
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.07699 -0.02305 -0.00018  0.02396  0.09492 
## 
## Coefficients:
##                                  Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                     1.793e-02  1.137e-03  15.764   <2e-16 ***
## pct_esp_vert_diff_2011.2019b500 2.744e-05  1.696e-04   0.162    0.872    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.02906 on 948 degrees of freedom
##   (20 observations deleted due to missingness)
## Multiple R-squared:  2.761e-05,  Adjusted R-squared:  -0.001027 
## F-statistic: 0.02618 on 1 and 948 DF,  p-value: 0.8715

2.3.1.4 Buffers 750m

2.3.1.4.1 Pampalon / Material vs. bike lane length change

ggplot(units::drop_units(bei_df), aes(x=Bike_lane_diff.2011.2016b750, y=wSCOREMAT.2016)) +
  geom_point() +
  geom_smooth(method=lm)

res.lm <- lm(wSCOREMAT.2016 ~ Bike_lane_diff.2011.2016b750, data = units::drop_units(bei_df))
summary(res.lm)

## 
## Call:
## lm(formula = wSCOREMAT.2016 ~ Bike_lane_diff.2011.2016b750, data = units::drop_units(bei_df))
## 
## Residuals:
##       Min        1Q    Median        3Q       Max 
## -0.142411 -0.022130 -0.000503  0.022855  0.137007 
## 
## Coefficients:
##                                Estimate Std. Error t value Pr(>|t|)   
## (Intercept)                   1.869e-03  1.372e-03   1.362  0.17341   
## Bike_lane_diff.2011.2016b750 -1.258e-06  3.904e-07  -3.221  0.00132 **
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.03594 on 948 degrees of freedom
##   (20 observations deleted due to missingness)
## Multiple R-squared:  0.01083,    Adjusted R-squared:  0.009782 
## F-statistic: 10.38 on 1 and 948 DF,  p-value: 0.001321

2.3.1.4.2 Pampalon / Social vs. bike lane length change

ggplot(units::drop_units(bei_df), aes(x=Bike_lane_diff.2011.2016b750, y=wSCORESOC.2016)) +
  geom_point() +
  geom_smooth(method=lm)

res.lm <- lm(wSCORESOC.2016 ~ Bike_lane_diff.2011.2016b750, data = units::drop_units(bei_df))
summary(res.lm)

## 
## Call:
## lm(formula = wSCORESOC.2016 ~ Bike_lane_diff.2011.2016b750, data = units::drop_units(bei_df))
## 
## Residuals:
##       Min        1Q    Median        3Q       Max 
## -0.089263 -0.020055  0.000102  0.020760  0.099702 
## 
## Coefficients:
##                               Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                  1.327e-02  1.070e-03  12.398   <2e-16 ***
## Bike_lane_diff.2011.2016b750 2.573e-06  3.045e-07   8.449   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.02803 on 948 degrees of freedom
##   (20 observations deleted due to missingness)
## Multiple R-squared:  0.07003,    Adjusted R-squared:  0.06905 
## F-statistic: 71.38 on 1 and 948 DF,  p-value: < 2.2e-16

2.3.1.4.3 Pampalon / Material vs. green space change

ggplot(units::drop_units(bei_df), aes(x=pct_esp_vert_diff_2011.2019b750, y=wSCOREMAT.2016)) +
  geom_point() +
  geom_smooth(method=lm)

res.lm <- lm(wSCOREMAT.2016 ~ pct_esp_vert_diff_2011.2019b750, data = units::drop_units(bei_df))
summary(res.lm)

## 
## Call:
## lm(formula = wSCOREMAT.2016 ~ pct_esp_vert_diff_2011.2019b750, 
##     data = units::drop_units(bei_df))
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.13978 -0.02219  0.00027  0.02354  0.13981 
## 
## Coefficients:
##                                   Estimate Std. Error t value Pr(>|t|)  
## (Intercept)                      0.0011957  0.0014218   0.841   0.4006  
## pct_esp_vert_diff_2011.2019b750 -0.0004450  0.0002172  -2.049   0.0407 *
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.03605 on 948 degrees of freedom
##   (20 observations deleted due to missingness)
## Multiple R-squared:  0.004409,   Adjusted R-squared:  0.003359 
## F-statistic: 4.198 on 1 and 948 DF,  p-value: 0.04074

2.3.1.4.4 Pampalon / Social vs. green space change

ggplot(units::drop_units(bei_df), aes(x=pct_esp_vert_diff_2011.2019b750, y=wSCORESOC.2016)) +
  geom_point() +
  geom_smooth(method=lm)

res.lm <- lm(wSCORESOC.2016 ~ pct_esp_vert_diff_2011.2019b750, data = units::drop_units(bei_df))
summary(res.lm)

## 
## Call:
## lm(formula = wSCORESOC.2016 ~ pct_esp_vert_diff_2011.2019b750, 
##     data = units::drop_units(bei_df))
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.07653 -0.02310 -0.00010  0.02383  0.09496 
## 
## Coefficients:
##                                   Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                      1.817e-02  1.146e-03   15.85   <2e-16 ***
## pct_esp_vert_diff_2011.2019b750 -3.675e-05  1.751e-04   -0.21    0.834    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.02906 on 948 degrees of freedom
##   (20 observations deleted due to missingness)
## Multiple R-squared:  4.647e-05,  Adjusted R-squared:  -0.001008 
## F-statistic: 0.04405 on 1 and 948 DF,  p-value: 0.8338

2.3.2 Gentrified CT outcome

2.3.2.1 Census Tracts

2.3.2.1.1 Gentrified CT in 2016 vs. bike lane length change

ggplot(drop_na(units::drop_units(bei_df), gentrified_2016_2011), aes(x=Bike_lane_diff.2011.2016ct, y=gentrified_2016_2011)) +
  geom_boxplot()

bei_df %>% 
  units::drop_units() %>%
  drop_na() %>%
  group_by(gentrified_2016_2011) %>%
  summarise(
    count = n(),
    mean = mean(Bike_lane_diff.2011.2016ct, na.rm = TRUE),
    sd = sd(Bike_lane_diff.2011.2016ct, na.rm = TRUE)
  )

# Compute the analysis of variance
res.aov <- aov(Bike_lane_diff.2011.2016ct ~ gentrified_2016_2011, data = units::drop_units(bei_df))
# Summary of the analysis
summary(res.aov)

##                       Df    Sum Sq Mean Sq F value Pr(>F)
## gentrified_2016_2011   1    205247  205247   0.253  0.615
## Residuals            945 765953237  810533               
## 23 observations deleted due to missingness

2.3.2.1.2 Gentrified CT in 2011 vs. bike lane length change

ggplot(drop_na(units::drop_units(bei_df), gentrified_2011_2006), aes(x=Bike_lane_diff.2011.2016ct, y=gentrified_2011_2006)) +
  geom_boxplot() 

bei_df %>% 
  units::drop_units() %>%
  drop_na() %>%
  group_by(gentrified_2011_2006) %>%
  summarise(
    count = n(),
    mean = mean(Bike_lane_diff.2011.2016ct, na.rm = TRUE),
    sd = sd(Bike_lane_diff.2011.2016ct, na.rm = TRUE)
  )

# Compute the analysis of variance
res.aov <- aov(Bike_lane_diff.2011.2016ct ~ gentrified_2011_2006, data = units::drop_units(bei_df))
# Summary of the analysis
summary(res.aov)

##                       Df    Sum Sq Mean Sq F value Pr(>F)
## gentrified_2011_2006   1    691401  691401   0.805   0.37
## Residuals            886 760756987  858642               
## 82 observations deleted due to missingness

2.3.2.1.3 Gentrified CT in 2006 vs. bike lane length change

ggplot(drop_na(units::drop_units(bei_df), gentrified_2006_2001), aes(x=Bike_lane_diff.2011.2016ct, y=gentrified_2006_2001)) +
  geom_boxplot() 

bei_df %>% 
  units::drop_units() %>%
  drop_na() %>%
  group_by(gentrified_2006_2001) %>%
  summarise(
    count = n(),
    mean = mean(Bike_lane_diff.2011.2016ct, na.rm = TRUE),
    sd = sd(Bike_lane_diff.2011.2016ct, na.rm = TRUE)
  )

# Compute the analysis of variance
res.aov <- aov(Bike_lane_diff.2011.2016ct ~ gentrified_2006_2001, data = units::drop_units(bei_df))
# Summary of the analysis
summary(res.aov)

##                       Df    Sum Sq Mean Sq F value Pr(>F)  
## gentrified_2006_2001   1   1629828 1629828   3.598 0.0582 .
## Residuals            802 363318022  453015                 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 166 observations deleted due to missingness

2.3.2.1.4 Gentrified CT in 2016 vs. green space change

ggplot(drop_na(units::drop_units(bei_df), gentrified_2016_2011), aes(x=pct_esp_vert_diff_2011.2019ct, y=gentrified_2016_2011)) +
  geom_boxplot() 

bei_df %>% 
  units::drop_units() %>%
  drop_na() %>%
  group_by(gentrified_2016_2011) %>%
  summarise(
    count = n(),
    mean = mean(pct_esp_vert_diff_2011.2019ct, na.rm = TRUE),
    sd = sd(pct_esp_vert_diff_2011.2019ct, na.rm = TRUE)
  )

# Compute the analysis of variance
res.aov <- aov(pct_esp_vert_diff_2011.2019ct ~ gentrified_2016_2011, data = units::drop_units(bei_df))
# Summary of the analysis
summary(res.aov)

##                       Df Sum Sq Mean Sq F value Pr(>F)
## gentrified_2016_2011   1     42   42.47   1.218   0.27
## Residuals            945  32955   34.87               
## 23 observations deleted due to missingness

2.3.2.1.5 Gentrified CT in 2011 vs. green space change

ggplot(drop_na(units::drop_units(bei_df), gentrified_2011_2006), aes(x=pct_esp_vert_diff_2011.2019ct, y=gentrified_2011_2006)) +
  geom_boxplot() 

bei_df %>% 
  units::drop_units() %>%
  drop_na() %>%
  group_by(gentrified_2011_2006) %>%
  summarise(
    count = n(),
    mean = mean(pct_esp_vert_diff_2011.2019ct, na.rm = TRUE),
    sd = sd(pct_esp_vert_diff_2011.2019ct, na.rm = TRUE)
  )

# Compute the analysis of variance
res.aov <- aov(pct_esp_vert_diff_2011.2019ct ~ gentrified_2011_2006, data = units::drop_units(bei_df))
# Summary of the analysis
summary(res.aov)

##                       Df Sum Sq Mean Sq F value Pr(>F)  
## gentrified_2011_2006   1    165  164.99   4.791 0.0289 *
## Residuals            886  30509   34.43                 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 82 observations deleted due to missingness

2.3.2.1.6 Gentrified CT in 2006 vs. green space change

ggplot(drop_na(units::drop_units(bei_df), gentrified_2006_2001), aes(x=pct_esp_vert_diff_2011.2019ct, y=gentrified_2006_2001)) +
  geom_boxplot() 

bei_df %>% 
  units::drop_units() %>%
  drop_na() %>%
  group_by(gentrified_2006_2001) %>%
  summarise(
    count = n(),
    mean = mean(pct_esp_vert_diff_2011.2019ct, na.rm = TRUE),
    sd = sd(pct_esp_vert_diff_2011.2019ct, na.rm = TRUE)
  )

# Compute the analysis of variance
res.aov <- aov(pct_esp_vert_diff_2011.2019ct ~ gentrified_2006_2001, data = units::drop_units(bei_df))
# Summary of the analysis
summary(res.aov)

##                       Df Sum Sq Mean Sq F value Pr(>F)
## gentrified_2006_2001   1     81   81.36   2.674  0.102
## Residuals            802  24399   30.42               
## 166 observations deleted due to missingness

2.3.2.2 Buffers 250m

2.3.2.2.1 Gentrified CT in 2016 vs. bike lane length change

ggplot(drop_na(units::drop_units(bei_df), gentrified_2016_2011), aes(x=Bike_lane_diff.2011.2016b250, y=gentrified_2016_2011)) +
  geom_boxplot()

bei_df %>% 
  units::drop_units() %>%
  drop_na() %>%
  group_by(gentrified_2016_2011) %>%
  summarise(
    count = n(),
    mean = mean(Bike_lane_diff.2011.2016b250, na.rm = TRUE),
    sd = sd(Bike_lane_diff.2011.2016b250, na.rm = TRUE)
  )

# Compute the analysis of variance
res.aov <- aov(Bike_lane_diff.2011.2016b250 ~ gentrified_2016_2011, data = units::drop_units(bei_df))
# Summary of the analysis
summary(res.aov)

##                       Df    Sum Sq  Mean Sq F value   Pr(>F)    
## gentrified_2016_2011   1 2.977e+07 29772990   13.79 0.000216 ***
## Residuals            945 2.040e+09  2159034                     
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 23 observations deleted due to missingness

2.3.2.2.2 Gentrified CT in 2011 vs. bike lane length change

ggplot(drop_na(units::drop_units(bei_df), gentrified_2011_2006), aes(x=Bike_lane_diff.2011.2016b250, y=gentrified_2011_2006)) +
  geom_boxplot() 

bei_df %>% 
  units::drop_units() %>%
  drop_na() %>%
  group_by(gentrified_2011_2006) %>%
  summarise(
    count = n(),
    mean = mean(Bike_lane_diff.2011.2016b250, na.rm = TRUE),
    sd = sd(Bike_lane_diff.2011.2016b250, na.rm = TRUE)
  )

# Compute the analysis of variance
res.aov <- aov(Bike_lane_diff.2011.2016b250 ~ gentrified_2011_2006, data = units::drop_units(bei_df))
# Summary of the analysis
summary(res.aov)

##                       Df    Sum Sq  Mean Sq F value   Pr(>F)    
## gentrified_2011_2006   1 3.867e+07 38666795   17.08 3.92e-05 ***
## Residuals            886 2.005e+09  2263514                     
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 82 observations deleted due to missingness

2.3.2.2.3 Gentrified CT in 2006 vs. bike lane length change

ggplot(drop_na(units::drop_units(bei_df), gentrified_2006_2001), aes(x=Bike_lane_diff.2011.2016b250, y=gentrified_2006_2001)) +
  geom_boxplot() 

bei_df %>% 
  units::drop_units() %>%
  drop_na() %>%
  group_by(gentrified_2006_2001) %>%
  summarise(
    count = n(),
    mean = mean(Bike_lane_diff.2011.2016b250, na.rm = TRUE),
    sd = sd(Bike_lane_diff.2011.2016b250, na.rm = TRUE)
  )

# Compute the analysis of variance
res.aov <- aov(Bike_lane_diff.2011.2016b250 ~ gentrified_2006_2001, data = units::drop_units(bei_df))
# Summary of the analysis
summary(res.aov)

##                       Df   Sum Sq  Mean Sq F value   Pr(>F)    
## gentrified_2006_2001   1 5.52e+07 55201554   40.24 3.74e-10 ***
## Residuals            802 1.10e+09  1371697                     
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 166 observations deleted due to missingness

2.3.2.2.4 Gentrified CT in 2016 vs. green space change

ggplot(drop_na(units::drop_units(bei_df), gentrified_2016_2011), aes(x=pct_esp_vert_diff_2011.2019b250, y=gentrified_2016_2011)) +
  geom_boxplot() 

bei_df %>% 
  units::drop_units() %>%
  drop_na() %>%
  group_by(gentrified_2016_2011) %>%
  summarise(
    count = n(),
    mean = mean(pct_esp_vert_diff_2011.2019b250, na.rm = TRUE),
    sd = sd(pct_esp_vert_diff_2011.2019b250, na.rm = TRUE)
  )

# Compute the analysis of variance
res.aov <- aov(pct_esp_vert_diff_2011.2019b250 ~ gentrified_2016_2011, data = units::drop_units(bei_df))
# Summary of the analysis
summary(res.aov)

##                       Df Sum Sq Mean Sq F value Pr(>F)
## gentrified_2016_2011   1     17   17.34    0.54  0.463
## Residuals            945  30337   32.10               
## 23 observations deleted due to missingness

2.3.2.2.5 Gentrified CT in 2011 vs. green space change

ggplot(drop_na(units::drop_units(bei_df), gentrified_2011_2006), aes(x=pct_esp_vert_diff_2011.2019b250, y=gentrified_2011_2006)) +
  geom_boxplot() 

bei_df %>% 
  units::drop_units() %>%
  drop_na() %>%
  group_by(gentrified_2011_2006) %>%
  summarise(
    count = n(),
    mean = mean(pct_esp_vert_diff_2011.2019b250, na.rm = TRUE),
    sd = sd(pct_esp_vert_diff_2011.2019b250, na.rm = TRUE)
  )

# Compute the analysis of variance
res.aov <- aov(pct_esp_vert_diff_2011.2019b250 ~ gentrified_2011_2006, data = units::drop_units(bei_df))
# Summary of the analysis
summary(res.aov)

##                       Df Sum Sq Mean Sq F value Pr(>F)  
## gentrified_2011_2006   1    147   147.0   5.367 0.0208 *
## Residuals            886  24275    27.4                 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 82 observations deleted due to missingness

2.3.2.2.6 Gentrified CT in 2006 vs. green space change

ggplot(drop_na(units::drop_units(bei_df), gentrified_2006_2001), aes(x=pct_esp_vert_diff_2011.2019b250, y=gentrified_2006_2001)) +
  geom_boxplot() 

bei_df %>% 
  units::drop_units() %>%
  drop_na() %>%
  group_by(gentrified_2006_2001) %>%
  summarise(
    count = n(),
    mean = mean(pct_esp_vert_diff_2011.2019b250, na.rm = TRUE),
    sd = sd(pct_esp_vert_diff_2011.2019b250, na.rm = TRUE)
  )

# Compute the analysis of variance
res.aov <- aov(pct_esp_vert_diff_2011.2019b250 ~ gentrified_2006_2001, data = units::drop_units(bei_df))
# Summary of the analysis
summary(res.aov)

##                       Df Sum Sq Mean Sq F value Pr(>F)  
## gentrified_2006_2001   1     87   86.81   3.495 0.0619 .
## Residuals            802  19923   24.84                 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 166 observations deleted due to missingness

2.3.2.3 Buffers 500m

2.3.2.3.1 Gentrified CT in 2016 vs. bike lane length change

ggplot(drop_na(units::drop_units(bei_df), gentrified_2016_2011), aes(x=Bike_lane_diff.2011.2016b500, y=gentrified_2016_2011)) +
  geom_boxplot()

bei_df %>% 
  units::drop_units() %>%
  drop_na() %>%
  group_by(gentrified_2016_2011) %>%
  summarise(
    count = n(),
    mean = mean(Bike_lane_diff.2011.2016b500, na.rm = TRUE),
    sd = sd(Bike_lane_diff.2011.2016b500, na.rm = TRUE)
  )

# Compute the analysis of variance
res.aov <- aov(Bike_lane_diff.2011.2016b500 ~ gentrified_2016_2011, data = units::drop_units(bei_df))
# Summary of the analysis
summary(res.aov)

##                       Df    Sum Sq   Mean Sq F value   Pr(>F)    
## gentrified_2016_2011   1 1.594e+08 159408890   35.28 4.01e-09 ***
## Residuals            945 4.270e+09   4518550                     
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 23 observations deleted due to missingness

2.3.2.3.2 Gentrified CT in 2011 vs. bike lane length change

ggplot(drop_na(units::drop_units(bei_df), gentrified_2011_2006), aes(x=Bike_lane_diff.2011.2016b500, y=gentrified_2011_2006)) +
  geom_boxplot() 

bei_df %>% 
  units::drop_units() %>%
  drop_na() %>%
  group_by(gentrified_2011_2006) %>%
  summarise(
    count = n(),
    mean = mean(Bike_lane_diff.2011.2016b500, na.rm = TRUE),
    sd = sd(Bike_lane_diff.2011.2016b500, na.rm = TRUE)
  )

# Compute the analysis of variance
res.aov <- aov(Bike_lane_diff.2011.2016b500 ~ gentrified_2011_2006, data = units::drop_units(bei_df))
# Summary of the analysis
summary(res.aov)

##                       Df    Sum Sq   Mean Sq F value  Pr(>F)    
## gentrified_2011_2006   1 1.625e+08 162483870   34.41 6.3e-09 ***
## Residuals            886 4.184e+09   4722050                    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 82 observations deleted due to missingness

2.3.2.3.3 Gentrified CT in 2006 vs. bike lane length change

ggplot(drop_na(units::drop_units(bei_df), gentrified_2006_2001), aes(x=Bike_lane_diff.2011.2016b500, y=gentrified_2006_2001)) +
  geom_boxplot() 

bei_df %>% 
  units::drop_units() %>%
  drop_na() %>%
  group_by(gentrified_2006_2001) %>%
  summarise(
    count = n(),
    mean = mean(Bike_lane_diff.2011.2016b500, na.rm = TRUE),
    sd = sd(Bike_lane_diff.2011.2016b500, na.rm = TRUE)
  )

# Compute the analysis of variance
res.aov <- aov(Bike_lane_diff.2011.2016b500 ~ gentrified_2006_2001, data = units::drop_units(bei_df))
# Summary of the analysis
summary(res.aov)

##                       Df    Sum Sq   Mean Sq F value Pr(>F)    
## gentrified_2006_2001   1 2.423e+08 242253012   71.08 <2e-16 ***
## Residuals            802 2.733e+09   3408209                   
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 166 observations deleted due to missingness

2.3.2.3.4 Gentrified CT in 2016 vs. green space change

ggplot(drop_na(units::drop_units(bei_df), gentrified_2016_2011), aes(x=pct_esp_vert_diff_2011.2019b500, y=gentrified_2016_2011)) +
  geom_boxplot() 

bei_df %>% 
  units::drop_units() %>%
  drop_na() %>%
  group_by(gentrified_2016_2011) %>%
  summarise(
    count = n(),
    mean = mean(pct_esp_vert_diff_2011.2019b500, na.rm = TRUE),
    sd = sd(pct_esp_vert_diff_2011.2019b500, na.rm = TRUE)
  )

# Compute the analysis of variance
res.aov <- aov(pct_esp_vert_diff_2011.2019b500 ~ gentrified_2016_2011, data = units::drop_units(bei_df))
# Summary of the analysis
summary(res.aov)

##                       Df Sum Sq Mean Sq F value Pr(>F)
## gentrified_2016_2011   1      9   8.524   0.282  0.596
## Residuals            945  28608  30.273               
## 23 observations deleted due to missingness

2.3.2.3.5 Gentrified CT in 2011 vs. green space change

ggplot(drop_na(units::drop_units(bei_df), gentrified_2011_2006), aes(x=pct_esp_vert_diff_2011.2019b500, y=gentrified_2011_2006)) +
  geom_boxplot() 

bei_df %>% 
  units::drop_units() %>%
  drop_na() %>%
  group_by(gentrified_2011_2006) %>%
  summarise(
    count = n(),
    mean = mean(pct_esp_vert_diff_2011.2019b500, na.rm = TRUE),
    sd = sd(pct_esp_vert_diff_2011.2019b500, na.rm = TRUE)
  )

# Compute the analysis of variance
res.aov <- aov(pct_esp_vert_diff_2011.2019b500 ~ gentrified_2011_2006, data = units::drop_units(bei_df))
# Summary of the analysis
summary(res.aov)

##                       Df Sum Sq Mean Sq F value Pr(>F)  
## gentrified_2011_2006   1     96   95.68    3.71 0.0544 .
## Residuals            886  22853   25.79                 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 82 observations deleted due to missingness

2.3.2.3.6 Gentrified CT in 2006 vs. green space change

ggplot(drop_na(units::drop_units(bei_df), gentrified_2006_2001), aes(x=pct_esp_vert_diff_2011.2019b500, y=gentrified_2006_2001)) +
  geom_boxplot() 

bei_df %>% 
  units::drop_units() %>%
  drop_na() %>%
  group_by(gentrified_2006_2001) %>%
  summarise(
    count = n(),
    mean = mean(pct_esp_vert_diff_2011.2019b500, na.rm = TRUE),
    sd = sd(pct_esp_vert_diff_2011.2019b500, na.rm = TRUE)
  )

# Compute the analysis of variance
res.aov <- aov(pct_esp_vert_diff_2011.2019b500 ~ gentrified_2006_2001, data = units::drop_units(bei_df))
# Summary of the analysis
summary(res.aov)

##                       Df Sum Sq Mean Sq F value Pr(>F)
## gentrified_2006_2001   1     60   59.77   2.554   0.11
## Residuals            802  18768   23.40               
## 166 observations deleted due to missingness

2.3.2.4 Buffers 750m

2.3.2.4.1 Gentrified CT in 2016 vs. bike lane length change

ggplot(drop_na(units::drop_units(bei_df), gentrified_2016_2011), aes(x=Bike_lane_diff.2011.2016b750, y=gentrified_2016_2011)) +
  geom_boxplot()

bei_df %>% 
  units::drop_units() %>%
  drop_na() %>%
  group_by(gentrified_2016_2011) %>%
  summarise(
    count = n(),
    mean = mean(Bike_lane_diff.2011.2016b750, na.rm = TRUE),
    sd = sd(Bike_lane_diff.2011.2016b750, na.rm = TRUE)
  )

# Compute the analysis of variance
res.aov <- aov(Bike_lane_diff.2011.2016b750 ~ gentrified_2016_2011, data = units::drop_units(bei_df))
# Summary of the analysis
summary(res.aov)

##                       Df    Sum Sq   Mean Sq F value   Pr(>F)    
## gentrified_2016_2011   1 4.517e+08 451658967   53.32 6.02e-13 ***
## Residuals            945 8.005e+09   8471075                     
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 23 observations deleted due to missingness

2.3.2.4.2 Gentrified CT in 2011 vs. bike lane length change

ggplot(drop_na(units::drop_units(bei_df), gentrified_2011_2006), aes(x=Bike_lane_diff.2011.2016b750, y=gentrified_2011_2006)) +
  geom_boxplot() 

bei_df %>% 
  units::drop_units() %>%
  drop_na() %>%
  group_by(gentrified_2011_2006) %>%
  summarise(
    count = n(),
    mean = mean(Bike_lane_diff.2011.2016b750, na.rm = TRUE),
    sd = sd(Bike_lane_diff.2011.2016b750, na.rm = TRUE)
  )

# Compute the analysis of variance
res.aov <- aov(Bike_lane_diff.2011.2016b750 ~ gentrified_2011_2006, data = units::drop_units(bei_df))
# Summary of the analysis
summary(res.aov)

##                       Df    Sum Sq   Mean Sq F value   Pr(>F)    
## gentrified_2011_2006   1 5.376e+08 537556297   61.72 1.14e-14 ***
## Residuals            886 7.717e+09   8709779                     
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 82 observations deleted due to missingness

2.3.2.4.3 Gentrified CT in 2006 vs. bike lane length change

ggplot(drop_na(units::drop_units(bei_df), gentrified_2006_2001), aes(x=Bike_lane_diff.2011.2016b750, y=gentrified_2006_2001)) +
  geom_boxplot() 

bei_df %>% 
  units::drop_units() %>%
  drop_na() %>%
  group_by(gentrified_2006_2001) %>%
  summarise(
    count = n(),
    mean = mean(Bike_lane_diff.2011.2016b750, na.rm = TRUE),
    sd = sd(Bike_lane_diff.2011.2016b750, na.rm = TRUE)
  )

# Compute the analysis of variance
res.aov <- aov(Bike_lane_diff.2011.2016b750 ~ gentrified_2006_2001, data = units::drop_units(bei_df))
# Summary of the analysis
summary(res.aov)

##                       Df    Sum Sq   Mean Sq F value Pr(>F)    
## gentrified_2006_2001   1 7.300e+08 729955328   105.4 <2e-16 ***
## Residuals            802 5.555e+09   6926433                   
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 166 observations deleted due to missingness

2.3.2.4.4 Gentrified CT in 2016 vs. green space change

ggplot(drop_na(units::drop_units(bei_df), gentrified_2016_2011), aes(x=pct_esp_vert_diff_2011.2019b750, y=gentrified_2016_2011)) +
  geom_boxplot() 

bei_df %>% 
  units::drop_units() %>%
  drop_na() %>%
  group_by(gentrified_2016_2011) %>%
  summarise(
    count = n(),
    mean = mean(pct_esp_vert_diff_2011.2019b750, na.rm = TRUE),
    sd = sd(pct_esp_vert_diff_2011.2019b750, na.rm = TRUE)
  )

# Compute the analysis of variance
res.aov <- aov(pct_esp_vert_diff_2011.2019b750 ~ gentrified_2016_2011, data = units::drop_units(bei_df))
# Summary of the analysis
summary(res.aov)

##                       Df Sum Sq Mean Sq F value Pr(>F)
## gentrified_2016_2011   1      1   1.054   0.037  0.847
## Residuals            945  26789  28.349               
## 23 observations deleted due to missingness

2.3.2.4.5 Gentrified CT in 2011 vs. green space change

ggplot(drop_na(units::drop_units(bei_df), gentrified_2011_2006), aes(x=pct_esp_vert_diff_2011.2019b750, y=gentrified_2011_2006)) +
  geom_boxplot() 

bei_df %>% 
  units::drop_units() %>%
  drop_na() %>%
  group_by(gentrified_2011_2006) %>%
  summarise(
    count = n(),
    mean = mean(pct_esp_vert_diff_2011.2019b750, na.rm = TRUE),
    sd = sd(pct_esp_vert_diff_2011.2019b750, na.rm = TRUE)
  )

# Compute the analysis of variance
res.aov <- aov(pct_esp_vert_diff_2011.2019b750 ~ gentrified_2011_2006, data = units::drop_units(bei_df))
# Summary of the analysis
summary(res.aov)

##                       Df Sum Sq Mean Sq F value Pr(>F)
## gentrified_2011_2006   1     58   58.42   2.378  0.123
## Residuals            886  21768   24.57               
## 82 observations deleted due to missingness

2.3.2.4.6 Gentrified CT in 2006 vs. green space change

ggplot(drop_na(units::drop_units(bei_df), gentrified_2006_2001), aes(x=pct_esp_vert_diff_2011.2019b750, y=gentrified_2006_2001)) +
  geom_boxplot() 

bei_df %>% 
  units::drop_units() %>%
  drop_na() %>%
  group_by(gentrified_2006_2001) %>%
  summarise(
    count = n(),
    mean = mean(pct_esp_vert_diff_2011.2019b750, na.rm = TRUE),
    sd = sd(pct_esp_vert_diff_2011.2019b750, na.rm = TRUE)
  )

# Compute the analysis of variance
res.aov <- aov(pct_esp_vert_diff_2011.2019b750 ~ gentrified_2006_2001, data = units::drop_units(bei_df))
# Summary of the analysis
summary(res.aov)

##                       Df Sum Sq Mean Sq F value Pr(>F)
## gentrified_2006_2001   1     28   28.27   1.254  0.263
## Residuals            802  18076   22.54               
## 166 observations deleted due to missingness

2.4 Association between outcomes and independant variables | INTERACT study area

2.4.1 Pampalon outcome

2.4.1.1 Census Tracts

2.4.1.1.1 Pampalon / Material vs. bike lane length change

ggplot(units::drop_units(filter(bei_df, interact_aoi)), aes(x=Bike_lane_diff.2011.2016ct, y=wSCOREMAT.2016)) +
  geom_point() +
  geom_smooth(method=lm)

res.lm <- lm(wSCOREMAT.2016 ~ Bike_lane_diff.2011.2016ct, data = units::drop_units(filter(bei_df, interact_aoi)))
summary(res.lm)

## 
## Call:
## lm(formula = wSCOREMAT.2016 ~ Bike_lane_diff.2011.2016ct, data = units::drop_units(filter(bei_df, 
##     interact_aoi)))
## 
## Residuals:
##       Min        1Q    Median        3Q       Max 
## -0.142085 -0.024870  0.001266  0.026881  0.137334 
## 
## Coefficients:
##                              Estimate Std. Error t value Pr(>|t|)  
## (Intercept)                 1.543e-03  1.605e-03   0.961   0.3368  
## Bike_lane_diff.2011.2016ct -3.559e-06  1.462e-06  -2.434   0.0152 *
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.03979 on 689 degrees of freedom
##   (14 observations deleted due to missingness)
## Multiple R-squared:  0.008529,   Adjusted R-squared:  0.00709 
## F-statistic: 5.927 on 1 and 689 DF,  p-value: 0.01517

2.4.1.1.2 Pampalon / Social vs. bike lane length change

ggplot(units::drop_units(filter(bei_df, interact_aoi)), aes(x=Bike_lane_diff.2011.2016ct, y=wSCORESOC.2016)) +
  geom_point() +
  geom_smooth(method=lm)

res.lm <- lm(wSCORESOC.2016 ~ Bike_lane_diff.2011.2016ct, data = units::drop_units(filter(bei_df, interact_aoi)))
summary(res.lm)

## 
## Call:
## lm(formula = wSCORESOC.2016 ~ Bike_lane_diff.2011.2016ct, data = units::drop_units(filter(bei_df, 
##     interact_aoi)))
## 
## Residuals:
##       Min        1Q    Median        3Q       Max 
## -0.077632 -0.020726  0.002851  0.023379  0.088820 
## 
## Coefficients:
##                              Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                 2.415e-02  1.173e-03  20.593   <2e-16 ***
## Bike_lane_diff.2011.2016ct -2.085e-06  1.068e-06  -1.953   0.0513 .  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.02907 on 689 degrees of freedom
##   (14 observations deleted due to missingness)
## Multiple R-squared:  0.005504,   Adjusted R-squared:  0.00406 
## F-statistic: 3.813 on 1 and 689 DF,  p-value: 0.05126

2.4.1.1.3 Pampalon / Material vs. green space change

ggplot(units::drop_units(filter(bei_df, interact_aoi)), aes(x=pct_esp_vert_diff_2011.2019ct, y=wSCOREMAT.2016)) +
  geom_point() +
  geom_smooth(method=lm)

res.lm <- lm(wSCOREMAT.2016 ~ pct_esp_vert_diff_2011.2019ct, data = units::drop_units(filter(bei_df, interact_aoi)))
summary(res.lm)

## 
## Call:
## lm(formula = wSCOREMAT.2016 ~ pct_esp_vert_diff_2011.2019ct, 
##     data = units::drop_units(filter(bei_df, interact_aoi)))
## 
## Residuals:
##       Min        1Q    Median        3Q       Max 
## -0.140225 -0.025875  0.001734  0.027493  0.133529 
## 
## Coefficients:
##                                 Estimate Std. Error t value Pr(>|t|)
## (Intercept)                   -0.0011341  0.0019505  -0.581    0.561
## pct_esp_vert_diff_2011.2019ct  0.0003967  0.0003530   1.124    0.261
## 
## Residual standard error: 0.03992 on 689 degrees of freedom
##   (14 observations deleted due to missingness)
## Multiple R-squared:  0.00183,    Adjusted R-squared:  0.0003809 
## F-statistic: 1.263 on 1 and 689 DF,  p-value: 0.2615

2.4.1.1.4 Pampalon / Social vs. green space change

ggplot(units::drop_units(filter(bei_df, interact_aoi)), aes(x=pct_esp_vert_diff_2011.2019ct, y=wSCORESOC.2016)) +
  geom_point() +
  geom_smooth(method=lm)

res.lm <- lm(wSCORESOC.2016 ~ pct_esp_vert_diff_2011.2019ct, data = units::drop_units(filter(bei_df, interact_aoi)))
summary(res.lm)

## 
## Call:
## lm(formula = wSCORESOC.2016 ~ pct_esp_vert_diff_2011.2019ct, 
##     data = units::drop_units(filter(bei_df, interact_aoi)))
## 
## Residuals:
##       Min        1Q    Median        3Q       Max 
## -0.077431 -0.020217  0.003207  0.020227  0.091996 
## 
## Coefficients:
##                                Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                   0.0174407  0.0013774  12.662  < 2e-16 ***
## pct_esp_vert_diff_2011.2019ct 0.0017147  0.0002493   6.879 1.36e-11 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.02819 on 689 degrees of freedom
##   (14 observations deleted due to missingness)
## Multiple R-squared:  0.06426,    Adjusted R-squared:  0.0629 
## F-statistic: 47.32 on 1 and 689 DF,  p-value: 1.357e-11

2.4.1.2 Buffers 250m

2.4.1.2.1 Pampalon / Material vs. bike lane length change

ggplot(units::drop_units(filter(bei_df, interact_aoi)), aes(x=Bike_lane_diff.2011.2016b250, y=wSCOREMAT.2016)) +
  geom_point() +
  geom_smooth(method=lm)

res.lm <- lm(wSCOREMAT.2016 ~ Bike_lane_diff.2011.2016b250, data = units::drop_units(filter(bei_df, interact_aoi)))
summary(res.lm)

## 
## Call:
## lm(formula = wSCOREMAT.2016 ~ Bike_lane_diff.2011.2016b250, data = units::drop_units(filter(bei_df, 
##     interact_aoi)))
## 
## Residuals:
##       Min        1Q    Median        3Q       Max 
## -0.143165 -0.025609  0.002509  0.026293  0.136254 
## 
## Coefficients:
##                                Estimate Std. Error t value Pr(>|t|)   
## (Intercept)                   2.623e-03  1.716e-03   1.528  0.12694   
## Bike_lane_diff.2011.2016b250 -2.651e-06  9.069e-07  -2.924  0.00357 **
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.03972 on 689 degrees of freedom
##   (14 observations deleted due to missingness)
## Multiple R-squared:  0.01225,    Adjusted R-squared:  0.01082 
## F-statistic: 8.547 on 1 and 689 DF,  p-value: 0.003574

2.4.1.2.2 Pampalon / Social vs. bike lane length change

ggplot(units::drop_units(filter(bei_df, interact_aoi)), aes(x=Bike_lane_diff.2011.2016b250, y=wSCORESOC.2016)) +
  geom_point() +
  geom_smooth(method=lm)

res.lm <- lm(wSCORESOC.2016 ~ Bike_lane_diff.2011.2016b250, data = units::drop_units(filter(bei_df, interact_aoi)))
summary(res.lm)

## 
## Call:
## lm(formula = wSCORESOC.2016 ~ Bike_lane_diff.2011.2016b250, data = units::drop_units(filter(bei_df, 
##     interact_aoi)))
## 
## Residuals:
##       Min        1Q    Median        3Q       Max 
## -0.076323 -0.020636  0.003626  0.023255  0.090129 
## 
## Coefficients:
##                               Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                  2.284e-02  1.259e-03  18.142   <2e-16 ***
## Bike_lane_diff.2011.2016b250 6.078e-07  6.651e-07   0.914    0.361    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.02913 on 689 degrees of freedom
##   (14 observations deleted due to missingness)
## Multiple R-squared:  0.001211,   Adjusted R-squared:  -0.000239 
## F-statistic: 0.8351 on 1 and 689 DF,  p-value: 0.3611

2.4.1.2.3 Pampalon / Material vs. green space change

ggplot(units::drop_units(filter(bei_df, interact_aoi)), aes(x=pct_esp_vert_diff_2011.2019b250, y=wSCOREMAT.2016)) +
  geom_point() +
  geom_smooth(method=lm)

res.lm <- lm(wSCOREMAT.2016 ~ pct_esp_vert_diff_2011.2019b250, data = units::drop_units(filter(bei_df, interact_aoi)))
summary(res.lm)

## 
## Call:
## lm(formula = wSCOREMAT.2016 ~ pct_esp_vert_diff_2011.2019b250, 
##     data = units::drop_units(filter(bei_df, interact_aoi)))
## 
## Residuals:
##       Min        1Q    Median        3Q       Max 
## -0.140800 -0.026341  0.001769  0.027364  0.138323 
## 
## Coefficients:
##                                   Estimate Std. Error t value Pr(>|t|)
## (Intercept)                     -0.0001825  0.0019837  -0.092    0.927
## pct_esp_vert_diff_2011.2019b250  0.0001318  0.0003963   0.333    0.740
## 
## Residual standard error: 0.03996 on 689 degrees of freedom
##   (14 observations deleted due to missingness)
## Multiple R-squared:  0.0001604,  Adjusted R-squared:  -0.001291 
## F-statistic: 0.1106 on 1 and 689 DF,  p-value: 0.7396

2.4.1.2.4 Pampalon / Social vs. green space change

ggplot(units::drop_units(filter(bei_df, interact_aoi)), aes(x=pct_esp_vert_diff_2011.2019b250, y=wSCORESOC.2016)) +
  geom_point() +
  geom_smooth(method=lm)

res.lm <- lm(wSCORESOC.2016 ~ pct_esp_vert_diff_2011.2019b250, data = units::drop_units(filter(bei_df, interact_aoi)))
summary(res.lm)

## 
## Call:
## lm(formula = wSCORESOC.2016 ~ pct_esp_vert_diff_2011.2019b250, 
##     data = units::drop_units(filter(bei_df, interact_aoi)))
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.08941 -0.01937  0.00312  0.02027  0.08932 
## 
## Coefficients:
##                                  Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                     0.0164684  0.0013876  11.868  < 2e-16 ***
## pct_esp_vert_diff_2011.2019b250 0.0021508  0.0002772   7.759 3.09e-14 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.02795 on 689 degrees of freedom
##   (14 observations deleted due to missingness)
## Multiple R-squared:  0.08036,    Adjusted R-squared:  0.07902 
## F-statistic:  60.2 on 1 and 689 DF,  p-value: 3.092e-14

2.4.1.3 Buffers 500m

2.4.1.3.1 Pampalon / Material vs. bike lane length change

ggplot(units::drop_units(filter(bei_df, interact_aoi)), aes(x=Bike_lane_diff.2011.2016b500, y=wSCOREMAT.2016)) +
  geom_point() +
  geom_smooth(method=lm)

res.lm <- lm(wSCOREMAT.2016 ~ Bike_lane_diff.2011.2016b500, data = units::drop_units(filter(bei_df, interact_aoi)))
summary(res.lm)

## 
## Call:
## lm(formula = wSCOREMAT.2016 ~ Bike_lane_diff.2011.2016b500, data = units::drop_units(filter(bei_df, 
##     interact_aoi)))
## 
## Residuals:
##       Min        1Q    Median        3Q       Max 
## -0.144252 -0.025817  0.001971  0.026545  0.135167 
## 
## Coefficients:
##                                Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                   3.710e-03  1.819e-03   2.039 0.041793 *  
## Bike_lane_diff.2011.2016b500 -2.151e-06  6.312e-07  -3.407 0.000694 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.03963 on 689 degrees of freedom
##   (14 observations deleted due to missingness)
## Multiple R-squared:  0.01657,    Adjusted R-squared:  0.01515 
## F-statistic: 11.61 on 1 and 689 DF,  p-value: 0.0006938

2.4.1.3.2 Pampalon / Social vs. bike lane length change

ggplot(units::drop_units(filter(bei_df, interact_aoi)), aes(x=Bike_lane_diff.2011.2016b500, y=wSCORESOC.2016)) +
  geom_point() +
  geom_smooth(method=lm)

res.lm <- lm(wSCORESOC.2016 ~ Bike_lane_diff.2011.2016b500, data = units::drop_units(filter(bei_df, interact_aoi)))
summary(res.lm)

## 
## Call:
## lm(formula = wSCORESOC.2016 ~ Bike_lane_diff.2011.2016b500, data = units::drop_units(filter(bei_df, 
##     interact_aoi)))
## 
## Residuals:
##       Min        1Q    Median        3Q       Max 
## -0.075316 -0.020362  0.003956  0.022246  0.092087 
## 
## Coefficients:
##                               Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                  2.088e-02  1.327e-03  15.737  < 2e-16 ***
## Bike_lane_diff.2011.2016b500 1.553e-06  4.605e-07   3.372 0.000788 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.02891 on 689 degrees of freedom
##   (14 observations deleted due to missingness)
## Multiple R-squared:  0.01624,    Adjusted R-squared:  0.01481 
## F-statistic: 11.37 on 1 and 689 DF,  p-value: 0.0007877

2.4.1.3.3 Pampalon / Material vs. green space change

ggplot(units::drop_units(filter(bei_df, interact_aoi)), aes(x=pct_esp_vert_diff_2011.2019b500, y=wSCOREMAT.2016)) +
  geom_point() +
  geom_smooth(method=lm)

res.lm <- lm(wSCOREMAT.2016 ~ pct_esp_vert_diff_2011.2019b500, data = units::drop_units(filter(bei_df, interact_aoi)))
summary(res.lm)

## 
## Call:
## lm(formula = wSCOREMAT.2016 ~ pct_esp_vert_diff_2011.2019b500, 
##     data = units::drop_units(filter(bei_df, interact_aoi)))
## 
## Residuals:
##       Min        1Q    Median        3Q       Max 
## -0.141012 -0.026324  0.001867  0.027400  0.138192 
## 
## Coefficients:
##                                   Estimate Std. Error t value Pr(>|t|)
## (Intercept)                     -0.0003420  0.0020046  -0.171    0.865
## pct_esp_vert_diff_2011.2019b500  0.0001867  0.0004184   0.446    0.656
## 
## Residual standard error: 0.03995 on 689 degrees of freedom
##   (14 observations deleted due to missingness)
## Multiple R-squared:  0.0002889,  Adjusted R-squared:  -0.001162 
## F-statistic: 0.1991 on 1 and 689 DF,  p-value: 0.6556

2.4.1.3.4 Pampalon / Social vs. green space change

ggplot(units::drop_units(filter(bei_df, interact_aoi)), aes(x=pct_esp_vert_diff_2011.2019b500, y=wSCORESOC.2016)) +
  geom_point() +
  geom_smooth(method=lm)

res.lm <- lm(wSCORESOC.2016 ~ pct_esp_vert_diff_2011.2019b500, data = units::drop_units(filter(bei_df, interact_aoi)))
summary(res.lm)

## 
## Call:
## lm(formula = wSCORESOC.2016 ~ pct_esp_vert_diff_2011.2019b500, 
##     data = units::drop_units(filter(bei_df, interact_aoi)))
## 
## Residuals:
##       Min        1Q    Median        3Q       Max 
## -0.092709 -0.018960  0.003189  0.020466  0.086738 
## 
## Coefficients:
##                                  Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                     0.0161365  0.0013996  11.529  < 2e-16 ***
## pct_esp_vert_diff_2011.2019b500 0.0023205  0.0002921   7.944 7.98e-15 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.0279 on 689 degrees of freedom
##   (14 observations deleted due to missingness)
## Multiple R-squared:  0.08391,    Adjusted R-squared:  0.08258 
## F-statistic: 63.11 on 1 and 689 DF,  p-value: 7.975e-15

2.4.1.4 Buffers 750m

2.4.1.4.1 Pampalon / Material vs. bike lane length change

ggplot(units::drop_units(filter(bei_df, interact_aoi)), aes(x=Bike_lane_diff.2011.2016b750, y=wSCOREMAT.2016)) +
  geom_point() +
  geom_smooth(method=lm)

res.lm <- lm(wSCOREMAT.2016 ~ Bike_lane_diff.2011.2016b750, data = units::drop_units(filter(bei_df, interact_aoi)))
summary(res.lm)

## 
## Call:
## lm(formula = wSCOREMAT.2016 ~ Bike_lane_diff.2011.2016b750, data = units::drop_units(filter(bei_df, 
##     interact_aoi)))
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.14482 -0.02576  0.00209  0.02703  0.13460 
## 
## Coefficients:
##                                Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                   4.274e-03  1.911e-03   2.236 0.025643 *  
## Bike_lane_diff.2011.2016b750 -1.594e-06  4.645e-07  -3.432 0.000634 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.03962 on 689 degrees of freedom
##   (14 observations deleted due to missingness)
## Multiple R-squared:  0.01681,    Adjusted R-squared:  0.01538 
## F-statistic: 11.78 on 1 and 689 DF,  p-value: 0.0006341

2.4.1.4.2 Pampalon / Social vs. bike lane length change

ggplot(units::drop_units(filter(bei_df, interact_aoi)), aes(x=Bike_lane_diff.2011.2016b750, y=wSCORESOC.2016)) +
  geom_point() +
  geom_smooth(method=lm)

res.lm <- lm(wSCORESOC.2016 ~ Bike_lane_diff.2011.2016b750, data = units::drop_units(filter(bei_df, interact_aoi)))
summary(res.lm)

## 
## Call:
## lm(formula = wSCORESOC.2016 ~ Bike_lane_diff.2011.2016b750, data = units::drop_units(filter(bei_df, 
##     interact_aoi)))
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.07989 -0.02017  0.00368  0.02173  0.09403 
## 
## Coefficients:
##                               Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                  1.894e-02  1.379e-03  13.738  < 2e-16 ***
## Bike_lane_diff.2011.2016b750 1.758e-06  3.350e-07   5.247 2.06e-07 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.02858 on 689 degrees of freedom
##   (14 observations deleted due to missingness)
## Multiple R-squared:  0.03843,    Adjusted R-squared:  0.03703 
## F-statistic: 27.54 on 1 and 689 DF,  p-value: 2.055e-07

2.4.1.4.3 Pampalon / Material vs. green space change

ggplot(units::drop_units(filter(bei_df, interact_aoi)), aes(x=pct_esp_vert_diff_2011.2019b750, y=wSCOREMAT.2016)) +
  geom_point() +
  geom_smooth(method=lm)

res.lm <- lm(wSCOREMAT.2016 ~ pct_esp_vert_diff_2011.2019b750, data = units::drop_units(filter(bei_df, interact_aoi)))
summary(res.lm)

## 
## Call:
## lm(formula = wSCOREMAT.2016 ~ pct_esp_vert_diff_2011.2019b750, 
##     data = units::drop_units(filter(bei_df, interact_aoi)))
## 
## Residuals:
##       Min        1Q    Median        3Q       Max 
## -0.140898 -0.026431  0.001807  0.027269  0.138488 
## 
## Coefficients:
##                                   Estimate Std. Error t value Pr(>|t|)
## (Intercept)                     -1.795e-05  2.012e-03  -0.009    0.993
## pct_esp_vert_diff_2011.2019b750  8.494e-05  4.317e-04   0.197    0.844
## 
## Residual standard error: 0.03996 on 689 degrees of freedom
##   (14 observations deleted due to missingness)
## Multiple R-squared:  5.618e-05,  Adjusted R-squared:  -0.001395 
## F-statistic: 0.03871 on 1 and 689 DF,  p-value: 0.8441

2.4.1.4.4 Pampalon / Social vs. green space change

ggplot(units::drop_units(filter(bei_df, interact_aoi)), aes(x=pct_esp_vert_diff_2011.2019b750, y=wSCORESOC.2016)) +
  geom_point() +
  geom_smooth(method=lm)

res.lm <- lm(wSCORESOC.2016 ~ pct_esp_vert_diff_2011.2019b750, data = units::drop_units(filter(bei_df, interact_aoi)))
summary(res.lm)

## 
## Call:
## lm(formula = wSCORESOC.2016 ~ pct_esp_vert_diff_2011.2019b750, 
##     data = units::drop_units(filter(bei_df, interact_aoi)))
## 
## Residuals:
##       Min        1Q    Median        3Q       Max 
## -0.087277 -0.019561  0.003355  0.020851  0.086431 
## 
## Coefficients:
##                                 Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                     0.016253   0.001407  11.549  < 2e-16 ***
## pct_esp_vert_diff_2011.2019b750 0.002338   0.000302   7.739 3.57e-14 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.02796 on 689 degrees of freedom
##   (14 observations deleted due to missingness)
## Multiple R-squared:  0.07998,    Adjusted R-squared:  0.07865 
## F-statistic:  59.9 on 1 and 689 DF,  p-value: 3.566e-14

2.4.2 Gentrified CT outcome

2.4.2.1 Census Tracts

2.4.2.1.1 Gentrified CT in 2016 vs. bike lane length change

ggplot(drop_na(units::drop_units(filter(bei_df, interact_aoi)), gentrified_2016_2011), aes(x=Bike_lane_diff.2011.2016ct, y=gentrified_2016_2011)) +
  geom_boxplot()

filter(bei_df, interact_aoi) %>% 
  units::drop_units() %>%
  drop_na() %>%
  group_by(gentrified_2016_2011) %>%
  summarise(
    count = n(),
    mean = mean(Bike_lane_diff.2011.2016ct, na.rm = TRUE),
    sd = sd(Bike_lane_diff.2011.2016ct, na.rm = TRUE)
  )

# Compute the analysis of variance
res.aov <- aov(Bike_lane_diff.2011.2016ct ~ gentrified_2016_2011, data = units::drop_units(filter(bei_df, interact_aoi)))
# Summary of the analysis
summary(res.aov)

##                       Df    Sum Sq Mean Sq F value Pr(>F)
## gentrified_2016_2011   1   1115825 1115825   1.038  0.309
## Residuals            688 739651503 1075075               
## 15 observations deleted due to missingness

2.4.2.1.2 Gentrified CT in 2011 vs. bike lane length change

ggplot(drop_na(units::drop_units(filter(bei_df, interact_aoi)), gentrified_2011_2006), aes(x=Bike_lane_diff.2011.2016ct, y=gentrified_2011_2006)) +
  geom_boxplot() 

filter(bei_df, interact_aoi) %>% 
  units::drop_units() %>%
  drop_na() %>%
  group_by(gentrified_2011_2006) %>%
  summarise(
    count = n(),
    mean = mean(Bike_lane_diff.2011.2016ct, na.rm = TRUE),
    sd = sd(Bike_lane_diff.2011.2016ct, na.rm = TRUE)
  )

# Compute the analysis of variance
res.aov <- aov(Bike_lane_diff.2011.2016ct ~ gentrified_2011_2006, data = units::drop_units(filter(bei_df, interact_aoi)))
# Summary of the analysis
summary(res.aov)

##                       Df    Sum Sq Mean Sq F value Pr(>F)
## gentrified_2011_2006   1    227980  227980   0.208  0.648
## Residuals            674 738566819 1095796               
## 29 observations deleted due to missingness

2.4.2.1.3 Gentrified CT in 2006 vs. bike lane length change

ggplot(drop_na(units::drop_units(filter(bei_df, interact_aoi)), gentrified_2006_2001), aes(x=Bike_lane_diff.2011.2016ct, y=gentrified_2006_2001)) +
  geom_boxplot() 

filter(bei_df, interact_aoi) %>% 
  units::drop_units() %>%
  drop_na() %>%
  group_by(gentrified_2006_2001) %>%
  summarise(
    count = n(),
    mean = mean(Bike_lane_diff.2011.2016ct, na.rm = TRUE),
    sd = sd(Bike_lane_diff.2011.2016ct, na.rm = TRUE)
  )

# Compute the analysis of variance
res.aov <- aov(Bike_lane_diff.2011.2016ct ~ gentrified_2006_2001, data = units::drop_units(filter(bei_df, interact_aoi)))
# Summary of the analysis
summary(res.aov)

##                       Df    Sum Sq Mean Sq F value Pr(>F)
## gentrified_2006_2001   1     27515   27515    0.05  0.824
## Residuals            634 351091051  553771               
## 69 observations deleted due to missingness

2.4.2.1.4 Gentrified CT in 2016 vs. green space change

ggplot(drop_na(units::drop_units(filter(bei_df, interact_aoi)), gentrified_2016_2011), aes(x=pct_esp_vert_diff_2011.2019ct, y=gentrified_2016_2011)) +
  geom_boxplot() 

filter(bei_df, interact_aoi) %>% 
  units::drop_units() %>%
  drop_na() %>%
  group_by(gentrified_2016_2011) %>%
  summarise(
    count = n(),
    mean = mean(pct_esp_vert_diff_2011.2019ct, na.rm = TRUE),
    sd = sd(pct_esp_vert_diff_2011.2019ct, na.rm = TRUE)
  )

# Compute the analysis of variance
res.aov <- aov(pct_esp_vert_diff_2011.2019ct ~ gentrified_2016_2011, data = units::drop_units(filter(bei_df, interact_aoi)))
# Summary of the analysis
summary(res.aov)

##                       Df Sum Sq Mean Sq F value   Pr(>F)    
## gentrified_2016_2011   1    285  284.57   15.67 8.32e-05 ***
## Residuals            688  12494   18.16                     
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 15 observations deleted due to missingness

2.4.2.1.5 Gentrified CT in 2011 vs. green space change

ggplot(drop_na(units::drop_units(filter(bei_df, interact_aoi)), gentrified_2011_2006), aes(x=pct_esp_vert_diff_2011.2019ct, y=gentrified_2011_2006)) +
  geom_boxplot() 

filter(bei_df, interact_aoi) %>% 
  units::drop_units() %>%
  drop_na() %>%
  group_by(gentrified_2011_2006) %>%
  summarise(
    count = n(),
    mean = mean(pct_esp_vert_diff_2011.2019ct, na.rm = TRUE),
    sd = sd(pct_esp_vert_diff_2011.2019ct, na.rm = TRUE)
  )

# Compute the analysis of variance
res.aov <- aov(pct_esp_vert_diff_2011.2019ct ~ gentrified_2011_2006, data = units::drop_units(filter(bei_df, interact_aoi)))
# Summary of the analysis
summary(res.aov)

##                       Df Sum Sq Mean Sq F value   Pr(>F)    
## gentrified_2011_2006   1    474   474.1   26.19 4.04e-07 ***
## Residuals            674  12202    18.1                     
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 29 observations deleted due to missingness

2.4.2.1.6 Gentrified CT in 2006 vs. green space change

ggplot(drop_na(units::drop_units(filter(bei_df, interact_aoi)), gentrified_2006_2001), aes(x=pct_esp_vert_diff_2011.2019ct, y=gentrified_2006_2001)) +
  geom_boxplot() 

filter(bei_df, interact_aoi) %>% 
  units::drop_units() %>%
  drop_na() %>%
  group_by(gentrified_2006_2001) %>%
  summarise(
    count = n(),
    mean = mean(pct_esp_vert_diff_2011.2019ct, na.rm = TRUE),
    sd = sd(pct_esp_vert_diff_2011.2019ct, na.rm = TRUE)
  )

# Compute the analysis of variance
res.aov <- aov(pct_esp_vert_diff_2011.2019ct ~ gentrified_2006_2001, data = units::drop_units(filter(bei_df, interact_aoi)))
# Summary of the analysis
summary(res.aov)

##                       Df Sum Sq Mean Sq F value   Pr(>F)    
## gentrified_2006_2001   1    277  277.15   16.69 4.97e-05 ***
## Residuals            634  10530   16.61                     
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 69 observations deleted due to missingness

2.4.2.2 Buffers 250m

2.4.2.2.1 Gentrified CT in 2016 vs. bike lane length change

ggplot(drop_na(units::drop_units(filter(bei_df, interact_aoi)), gentrified_2016_2011), aes(x=Bike_lane_diff.2011.2016b250, y=gentrified_2016_2011)) +
  geom_boxplot()

filter(bei_df, interact_aoi) %>% 
  units::drop_units() %>%
  drop_na() %>%
  group_by(gentrified_2016_2011) %>%
  summarise(
    count = n(),
    mean = mean(Bike_lane_diff.2011.2016b250, na.rm = TRUE),
    sd = sd(Bike_lane_diff.2011.2016b250, na.rm = TRUE)
  )

# Compute the analysis of variance
res.aov <- aov(Bike_lane_diff.2011.2016b250 ~ gentrified_2016_2011, data = units::drop_units(filter(bei_df, interact_aoi)))
# Summary of the analysis
summary(res.aov)

##                       Df    Sum Sq Mean Sq F value Pr(>F)
## gentrified_2016_2011   1 4.340e+06 4340160    1.56  0.212
## Residuals            688 1.914e+09 2781403               
## 15 observations deleted due to missingness

2.4.2.2.2 Gentrified CT in 2011 vs. bike lane length change

ggplot(drop_na(units::drop_units(filter(bei_df, interact_aoi)), gentrified_2011_2006), aes(x=Bike_lane_diff.2011.2016b250, y=gentrified_2011_2006)) +
  geom_boxplot() 

filter(bei_df, interact_aoi) %>% 
  units::drop_units() %>%
  drop_na() %>%
  group_by(gentrified_2011_2006) %>%
  summarise(
    count = n(),
    mean = mean(Bike_lane_diff.2011.2016b250, na.rm = TRUE),
    sd = sd(Bike_lane_diff.2011.2016b250, na.rm = TRUE)
  )

# Compute the analysis of variance
res.aov <- aov(Bike_lane_diff.2011.2016b250 ~ gentrified_2011_2006, data = units::drop_units(filter(bei_df, interact_aoi)))
# Summary of the analysis
summary(res.aov)

##                       Df    Sum Sq  Mean Sq F value Pr(>F)  
## gentrified_2011_2006   1 1.097e+07 10973598   3.898 0.0487 *
## Residuals            674 1.897e+09  2815091                 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 29 observations deleted due to missingness

2.4.2.2.3 Gentrified CT in 2006 vs. bike lane length change

ggplot(drop_na(units::drop_units(filter(bei_df, interact_aoi)), gentrified_2006_2001), aes(x=Bike_lane_diff.2011.2016b250, y=gentrified_2006_2001)) +
  geom_boxplot() 

filter(bei_df, interact_aoi) %>% 
  units::drop_units() %>%
  drop_na() %>%
  group_by(gentrified_2006_2001) %>%
  summarise(
    count = n(),
    mean = mean(Bike_lane_diff.2011.2016b250, na.rm = TRUE),
    sd = sd(Bike_lane_diff.2011.2016b250, na.rm = TRUE)
  )

# Compute the analysis of variance
res.aov <- aov(Bike_lane_diff.2011.2016b250 ~ gentrified_2006_2001, data = units::drop_units(filter(bei_df, interact_aoi)))
# Summary of the analysis
summary(res.aov)

##                       Df    Sum Sq  Mean Sq F value  Pr(>F)    
## gentrified_2006_2001   1 2.192e+07 21920663   13.42 0.00027 ***
## Residuals            634 1.036e+09  1633304                    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 69 observations deleted due to missingness

2.4.2.2.4 Gentrified CT in 2016 vs. green space change

ggplot(drop_na(units::drop_units(filter(bei_df, interact_aoi)), gentrified_2016_2011), aes(x=pct_esp_vert_diff_2011.2019b250, y=gentrified_2016_2011)) +
  geom_boxplot() 

filter(bei_df, interact_aoi) %>% 
  units::drop_units() %>%
  drop_na() %>%
  group_by(gentrified_2016_2011) %>%
  summarise(
    count = n(),
    mean = mean(pct_esp_vert_diff_2011.2019b250, na.rm = TRUE),
    sd = sd(pct_esp_vert_diff_2011.2019b250, na.rm = TRUE)
  )

# Compute the analysis of variance
res.aov <- aov(pct_esp_vert_diff_2011.2019b250 ~ gentrified_2016_2011, data = units::drop_units(filter(bei_df, interact_aoi)))
# Summary of the analysis
summary(res.aov)

##                       Df Sum Sq Mean Sq F value   Pr(>F)    
## gentrified_2016_2011   1    193  192.85   13.32 0.000283 ***
## Residuals            688   9963   14.48                     
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 15 observations deleted due to missingness

2.4.2.2.5 Gentrified CT in 2011 vs. green space change

ggplot(drop_na(units::drop_units(filter(bei_df, interact_aoi)), gentrified_2011_2006), aes(x=pct_esp_vert_diff_2011.2019b250, y=gentrified_2011_2006)) +
  geom_boxplot() 

filter(bei_df, interact_aoi) %>% 
  units::drop_units() %>%
  drop_na() %>%
  group_by(gentrified_2011_2006) %>%
  summarise(
    count = n(),
    mean = mean(pct_esp_vert_diff_2011.2019b250, na.rm = TRUE),
    sd = sd(pct_esp_vert_diff_2011.2019b250, na.rm = TRUE)
  )

# Compute the analysis of variance
res.aov <- aov(pct_esp_vert_diff_2011.2019b250 ~ gentrified_2011_2006, data = units::drop_units(filter(bei_df, interact_aoi)))
# Summary of the analysis
summary(res.aov)

##                       Df Sum Sq Mean Sq F value   Pr(>F)    
## gentrified_2011_2006   1    418   417.5   29.07 9.64e-08 ***
## Residuals            674   9679    14.4                     
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 29 observations deleted due to missingness

2.4.2.2.6 Gentrified CT in 2006 vs. green space change

ggplot(drop_na(units::drop_units(filter(bei_df, interact_aoi)), gentrified_2006_2001), aes(x=pct_esp_vert_diff_2011.2019b250, y=gentrified_2006_2001)) +
  geom_boxplot() 

filter(bei_df, interact_aoi) %>% 
  units::drop_units() %>%
  drop_na() %>%
  group_by(gentrified_2006_2001) %>%
  summarise(
    count = n(),
    mean = mean(pct_esp_vert_diff_2011.2019b250, na.rm = TRUE),
    sd = sd(pct_esp_vert_diff_2011.2019b250, na.rm = TRUE)
  )

# Compute the analysis of variance
res.aov <- aov(pct_esp_vert_diff_2011.2019b250 ~ gentrified_2006_2001, data = units::drop_units(filter(bei_df, interact_aoi)))
# Summary of the analysis
summary(res.aov)

##                       Df Sum Sq Mean Sq F value   Pr(>F)    
## gentrified_2006_2001   1    281  280.89   21.37 4.59e-06 ***
## Residuals            634   8333   13.14                     
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 69 observations deleted due to missingness

2.4.2.3 Buffers 500m

2.4.2.3.1 Gentrified CT in 2016 vs. bike lane length change

ggplot(drop_na(units::drop_units(filter(bei_df, interact_aoi)), gentrified_2016_2011), aes(x=Bike_lane_diff.2011.2016b500, y=gentrified_2016_2011)) +
  geom_boxplot()

filter(bei_df, interact_aoi) %>% 
  units::drop_units() %>%
  drop_na() %>%
  group_by(gentrified_2016_2011) %>%
  summarise(
    count = n(),
    mean = mean(Bike_lane_diff.2011.2016b500, na.rm = TRUE),
    sd = sd(Bike_lane_diff.2011.2016b500, na.rm = TRUE)
  )

# Compute the analysis of variance
res.aov <- aov(Bike_lane_diff.2011.2016b500 ~ gentrified_2016_2011, data = units::drop_units(filter(bei_df, interact_aoi)))
# Summary of the analysis
summary(res.aov)

##                       Df    Sum Sq  Mean Sq F value  Pr(>F)   
## gentrified_2016_2011   1 4.732e+07 47317317   8.356 0.00396 **
## Residuals            688 3.896e+09  5662417                   
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 15 observations deleted due to missingness

2.4.2.3.2 Gentrified CT in 2011 vs. bike lane length change

ggplot(drop_na(units::drop_units(filter(bei_df, interact_aoi)), gentrified_2011_2006), aes(x=Bike_lane_diff.2011.2016b500, y=gentrified_2011_2006)) +
  geom_boxplot() 

filter(bei_df, interact_aoi) %>% 
  units::drop_units() %>%
  drop_na() %>%
  group_by(gentrified_2011_2006) %>%
  summarise(
    count = n(),
    mean = mean(Bike_lane_diff.2011.2016b500, na.rm = TRUE),
    sd = sd(Bike_lane_diff.2011.2016b500, na.rm = TRUE)
  )

# Compute the analysis of variance
res.aov <- aov(Bike_lane_diff.2011.2016b500 ~ gentrified_2011_2006, data = units::drop_units(filter(bei_df, interact_aoi)))
# Summary of the analysis
summary(res.aov)

##                       Df    Sum Sq  Mean Sq F value  Pr(>F)   
## gentrified_2011_2006   1 6.005e+07 60051961   10.51 0.00125 **
## Residuals            674 3.852e+09  5715401                   
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 29 observations deleted due to missingness

2.4.2.3.3 Gentrified CT in 2006 vs. bike lane length change

ggplot(drop_na(units::drop_units(filter(bei_df, interact_aoi)), gentrified_2006_2001), aes(x=Bike_lane_diff.2011.2016b500, y=gentrified_2006_2001)) +
  geom_boxplot() 

filter(bei_df, interact_aoi) %>% 
  units::drop_units() %>%
  drop_na() %>%
  group_by(gentrified_2006_2001) %>%
  summarise(
    count = n(),
    mean = mean(Bike_lane_diff.2011.2016b500, na.rm = TRUE),
    sd = sd(Bike_lane_diff.2011.2016b500, na.rm = TRUE)
  )

# Compute the analysis of variance
res.aov <- aov(Bike_lane_diff.2011.2016b500 ~ gentrified_2006_2001, data = units::drop_units(filter(bei_df, interact_aoi)))
# Summary of the analysis
summary(res.aov)

##                       Df    Sum Sq   Mean Sq F value   Pr(>F)    
## gentrified_2006_2001   1 1.139e+08 113920313   28.58 1.25e-07 ***
## Residuals            634 2.527e+09   3985549                     
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 69 observations deleted due to missingness

2.4.2.3.4 Gentrified CT in 2016 vs. green space change

ggplot(drop_na(units::drop_units(filter(bei_df, interact_aoi)), gentrified_2016_2011), aes(x=pct_esp_vert_diff_2011.2019b500, y=gentrified_2016_2011)) +
  geom_boxplot() 

filter(bei_df, interact_aoi) %>% 
  units::drop_units() %>%
  drop_na() %>%
  group_by(gentrified_2016_2011) %>%
  summarise(
    count = n(),
    mean = mean(pct_esp_vert_diff_2011.2019b500, na.rm = TRUE),
    sd = sd(pct_esp_vert_diff_2011.2019b500, na.rm = TRUE)
  )

# Compute the analysis of variance
res.aov <- aov(pct_esp_vert_diff_2011.2019b500 ~ gentrified_2016_2011, data = units::drop_units(filter(bei_df, interact_aoi)))
# Summary of the analysis
summary(res.aov)

##                       Df Sum Sq Mean Sq F value   Pr(>F)    
## gentrified_2016_2011   1    162  161.57   12.41 0.000454 ***
## Residuals            688   8955   13.02                     
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 15 observations deleted due to missingness

2.4.2.3.5 Gentrified CT in 2011 vs. green space change

ggplot(drop_na(units::drop_units(filter(bei_df, interact_aoi)), gentrified_2011_2006), aes(x=pct_esp_vert_diff_2011.2019b500, y=gentrified_2011_2006)) +
  geom_boxplot() 

filter(bei_df, interact_aoi) %>% 
  units::drop_units() %>%
  drop_na() %>%
  group_by(gentrified_2011_2006) %>%
  summarise(
    count = n(),
    mean = mean(pct_esp_vert_diff_2011.2019b500, na.rm = TRUE),
    sd = sd(pct_esp_vert_diff_2011.2019b500, na.rm = TRUE)
  )

# Compute the analysis of variance
res.aov <- aov(pct_esp_vert_diff_2011.2019b500 ~ gentrified_2011_2006, data = units::drop_units(filter(bei_df, interact_aoi)))
# Summary of the analysis
summary(res.aov)

##                       Df Sum Sq Mean Sq F value   Pr(>F)    
## gentrified_2011_2006   1    357   356.9   27.58 2.03e-07 ***
## Residuals            674   8722    12.9                     
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 29 observations deleted due to missingness

2.4.2.3.6 Gentrified CT in 2006 vs. green space change

ggplot(drop_na(units::drop_units(filter(bei_df, interact_aoi)), gentrified_2006_2001), aes(x=pct_esp_vert_diff_2011.2019b500, y=gentrified_2006_2001)) +
  geom_boxplot() 

filter(bei_df, interact_aoi) %>% 
  units::drop_units() %>%
  drop_na() %>%
  group_by(gentrified_2006_2001) %>%
  summarise(
    count = n(),
    mean = mean(pct_esp_vert_diff_2011.2019b500, na.rm = TRUE),
    sd = sd(pct_esp_vert_diff_2011.2019b500, na.rm = TRUE)
  )

# Compute the analysis of variance
res.aov <- aov(pct_esp_vert_diff_2011.2019b500 ~ gentrified_2006_2001, data = units::drop_units(filter(bei_df, interact_aoi)))
# Summary of the analysis
summary(res.aov)

##                       Df Sum Sq Mean Sq F value   Pr(>F)    
## gentrified_2006_2001   1    245  244.66   20.61 6.74e-06 ***
## Residuals            634   7526   11.87                     
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 69 observations deleted due to missingness

2.4.2.4 Buffers 750m

2.4.2.4.1 Gentrified CT in 2016 vs. bike lane length change

ggplot(drop_na(units::drop_units(filter(bei_df, interact_aoi)), gentrified_2016_2011), aes(x=Bike_lane_diff.2011.2016b750, y=gentrified_2016_2011)) +
  geom_boxplot()

filter(bei_df, interact_aoi) %>% 
  units::drop_units() %>%
  drop_na() %>%
  group_by(gentrified_2016_2011) %>%
  summarise(
    count = n(),
    mean = mean(Bike_lane_diff.2011.2016b750, na.rm = TRUE),
    sd = sd(Bike_lane_diff.2011.2016b750, na.rm = TRUE)
  )

# Compute the analysis of variance
res.aov <- aov(Bike_lane_diff.2011.2016b750 ~ gentrified_2016_2011, data = units::drop_units(filter(bei_df, interact_aoi)))
# Summary of the analysis
summary(res.aov)

##                       Df    Sum Sq   Mean Sq F value   Pr(>F)    
## gentrified_2016_2011   1 1.527e+08 152740935   14.75 0.000134 ***
## Residuals            688 7.125e+09  10356102                     
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 15 observations deleted due to missingness

2.4.2.4.2 Gentrified CT in 2011 vs. bike lane length change

ggplot(drop_na(units::drop_units(filter(bei_df, interact_aoi)), gentrified_2011_2006), aes(x=Bike_lane_diff.2011.2016b750, y=gentrified_2011_2006)) +
  geom_boxplot() 

filter(bei_df, interact_aoi) %>% 
  units::drop_units() %>%
  drop_na() %>%
  group_by(gentrified_2011_2006) %>%
  summarise(
    count = n(),
    mean = mean(Bike_lane_diff.2011.2016b750, na.rm = TRUE),
    sd = sd(Bike_lane_diff.2011.2016b750, na.rm = TRUE)
  )

# Compute the analysis of variance
res.aov <- aov(Bike_lane_diff.2011.2016b750 ~ gentrified_2011_2006, data = units::drop_units(filter(bei_df, interact_aoi)))
# Summary of the analysis
summary(res.aov)

##                       Df    Sum Sq   Mean Sq F value   Pr(>F)    
## gentrified_2011_2006   1 2.429e+08 242939736   23.52 1.54e-06 ***
## Residuals            674 6.962e+09  10329115                     
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 29 observations deleted due to missingness

2.4.2.4.3 Gentrified CT in 2006 vs. bike lane length change

ggplot(drop_na(units::drop_units(filter(bei_df, interact_aoi)), gentrified_2006_2001), aes(x=Bike_lane_diff.2011.2016b750, y=gentrified_2006_2001)) +
  geom_boxplot() 

filter(bei_df, interact_aoi) %>% 
  units::drop_units() %>%
  drop_na() %>%
  group_by(gentrified_2006_2001) %>%
  summarise(
    count = n(),
    mean = mean(Bike_lane_diff.2011.2016b750, na.rm = TRUE),
    sd = sd(Bike_lane_diff.2011.2016b750, na.rm = TRUE)
  )

# Compute the analysis of variance
res.aov <- aov(Bike_lane_diff.2011.2016b750 ~ gentrified_2006_2001, data = units::drop_units(filter(bei_df, interact_aoi)))
# Summary of the analysis
summary(res.aov)

##                       Df    Sum Sq   Mean Sq F value  Pr(>F)    
## gentrified_2006_2001   1 3.754e+08 375368237   46.86 1.8e-11 ***
## Residuals            634 5.079e+09   8010555                    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 69 observations deleted due to missingness

2.4.2.4.4 Gentrified CT in 2016 vs. green space change

ggplot(drop_na(units::drop_units(filter(bei_df, interact_aoi)), gentrified_2016_2011), aes(x=pct_esp_vert_diff_2011.2019b750, y=gentrified_2016_2011)) +
  geom_boxplot() 

filter(bei_df, interact_aoi) %>% 
  units::drop_units() %>%
  drop_na() %>%
  group_by(gentrified_2016_2011) %>%
  summarise(
    count = n(),
    mean = mean(pct_esp_vert_diff_2011.2019b750, na.rm = TRUE),
    sd = sd(pct_esp_vert_diff_2011.2019b750, na.rm = TRUE)
  )

# Compute the analysis of variance
res.aov <- aov(pct_esp_vert_diff_2011.2019b750 ~ gentrified_2016_2011, data = units::drop_units(filter(bei_df, interact_aoi)))
# Summary of the analysis
summary(res.aov)

##                       Df Sum Sq Mean Sq F value  Pr(>F)   
## gentrified_2016_2011   1    133  133.34   10.88 0.00102 **
## Residuals            688   8432   12.26                   
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 15 observations deleted due to missingness

2.4.2.4.5 Gentrified CT in 2011 vs. green space change

ggplot(drop_na(units::drop_units(filter(bei_df, interact_aoi)), gentrified_2011_2006), aes(x=pct_esp_vert_diff_2011.2019b750, y=gentrified_2011_2006)) +
  geom_boxplot() 

filter(bei_df, interact_aoi) %>% 
  units::drop_units() %>%
  drop_na() %>%
  group_by(gentrified_2011_2006) %>%
  summarise(
    count = n(),
    mean = mean(pct_esp_vert_diff_2011.2019b750, na.rm = TRUE),
    sd = sd(pct_esp_vert_diff_2011.2019b750, na.rm = TRUE)
  )

# Compute the analysis of variance
res.aov <- aov(pct_esp_vert_diff_2011.2019b750 ~ gentrified_2011_2006, data = units::drop_units(filter(bei_df, interact_aoi)))
# Summary of the analysis
summary(res.aov)

##                       Df Sum Sq Mean Sq F value   Pr(>F)    
## gentrified_2011_2006   1    314  314.16   25.82 4.86e-07 ***
## Residuals            674   8200   12.17                     
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 29 observations deleted due to missingness

2.4.2.4.6 Gentrified CT in 2006 vs. green space change

ggplot(drop_na(units::drop_units(filter(bei_df, interact_aoi)), gentrified_2006_2001), aes(x=pct_esp_vert_diff_2011.2019b750, y=gentrified_2006_2001)) +
  geom_boxplot() 

filter(bei_df, interact_aoi) %>% 
  units::drop_units() %>%
  drop_na() %>%
  group_by(gentrified_2006_2001) %>%
  summarise(
    count = n(),
    mean = mean(pct_esp_vert_diff_2011.2019b750, na.rm = TRUE),
    sd = sd(pct_esp_vert_diff_2011.2019b750, na.rm = TRUE)
  )

# Compute the analysis of variance
res.aov <- aov(pct_esp_vert_diff_2011.2019b750 ~ gentrified_2006_2001, data = units::drop_units(filter(bei_df, interact_aoi)))
# Summary of the analysis
summary(res.aov)

##                       Df Sum Sq Mean Sq F value   Pr(>F)    
## gentrified_2006_2001   1    197  196.68   17.47 3.33e-05 ***
## Residuals            634   7139   11.26                     
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 69 observations deleted due to missingness

2.5 Association between outcomes and independant variables | Only CT with changes

2.5.1 Pampalon outcome

2.5.1.1 Census Tracts

ggplot(filter(units::drop_units(bei_df), Bike_lane_diff.2011.2016ct != 0), aes(x=Bike_lane_diff.2011.2016ct, y=wSCOREMAT.2016)) +
  geom_point() +
  geom_smooth(method=lm) +
  labs(title = "Pampalon / Material vs. bike lane length change")

ggplot(filter(units::drop_units(bei_df), Bike_lane_diff.2011.2016ct != 0), aes(x=Bike_lane_diff.2011.2016ct, y=wSCORESOC.2016)) +
  geom_point() +
  geom_smooth(method=lm) +
  labs(title = "Pampalon / Social vs. bike lane length change")

ggplot(filter(units::drop_units(bei_df), Bike_lane_diff.2011.2016ct != 0), aes(x=pct_esp_vert_diff_2011.2019ct, y=wSCOREMAT.2016)) +
  geom_point() +
  geom_smooth(method=lm) +
  labs(title = "Pampalon / Material vs. green space change")

ggplot(filter(units::drop_units(bei_df), Bike_lane_diff.2011.2016ct != 0), aes(x=pct_esp_vert_diff_2011.2019ct, y=wSCORESOC.2016)) +
  geom_point() +
  geom_smooth(method=lm) +
  labs(title = "Pampalon / Social vs. green space change")

2.5.1.2 Buffers 250m

TODO

2.5.1.3 Buffers 500m

TODO

2.5.1.4 Buffers 750m

TODO

2.5.2 Gentrified CT outcome

2.5.2.1 Census Tracts

ggplot(drop_na(filter(units::drop_units(bei_df), Bike_lane_diff.2011.2016ct != 0), gentrified_2016_2011), aes(x=Bike_lane_diff.2011.2016ct, y=gentrified_2016_2011)) +
  geom_boxplot() +
  labs(title = "Gentrified CT (2011 to 2016) vs. bike lane length change")

ggplot(drop_na(filter(units::drop_units(bei_df), Bike_lane_diff.2011.2016ct != 0), gentrified_2011_2006), aes(x=Bike_lane_diff.2011.2016ct, y=gentrified_2011_2006)) +
  geom_boxplot() +
  labs(title = "Gentrified CT (2006 to 2011) vs. bike lane length change")

ggplot(drop_na(filter(units::drop_units(bei_df), Bike_lane_diff.2011.2016ct != 0), gentrified_2006_2001), aes(x=Bike_lane_diff.2011.2016ct, y=gentrified_2006_2001)) +
  geom_boxplot() +
  labs(title = "Gentrified CT (2001 to 2006) vs. bike lane length change")

ggplot(drop_na(filter(units::drop_units(bei_df), Bike_lane_diff.2011.2016ct != 0), gentrified_2016_2011), aes(x=pct_esp_vert_diff_2011.2019ct, y=gentrified_2016_2011)) +
  geom_boxplot() +
  labs(title = "Gentrified CT (2011 to 2016) vs. green space change")

ggplot(drop_na(filter(units::drop_units(bei_df), Bike_lane_diff.2011.2016ct != 0), gentrified_2011_2006), aes(x=pct_esp_vert_diff_2011.2019ct, y=gentrified_2011_2006)) +
  geom_boxplot() +
  labs(title = "Gentrified CT (2006 to 2011) vs. green space change")

ggplot(drop_na(filter(units::drop_units(bei_df), Bike_lane_diff.2011.2016ct != 0), gentrified_2006_2001), aes(x=pct_esp_vert_diff_2011.2019ct, y=gentrified_2006_2001)) +
  geom_boxplot() +
  labs(title = "Gentrified CT (2001 to 2006) vs. green space change")

2.5.2.2 Buffers 250m

TODO

2.5.2.3 Buffers 500m

TODO

2.5.2.4 Buffers 750m

TODO